1. Introduction
A list of symbols is presented in Appendix A.
The nature of subglacial hydrology and basal motion on ice masses underlain by soft sediments are central questions in ice dynamics (e.g., Clarke, Reference Clarke1987; Murray, Reference Murray1997; Tulaczyk and others, Reference Tulaczyk, Kamb and Engelhardt2000). However, despite abundant evidence for subglacial sediments beneath fast-moving outlet glaciers and ice streams draining the Greenland and Antarctic ice sheets (e.g., Alley and others, Reference Alley, Blankenship, Bentley and Rooney1986; Blankenship and others, Reference Blankenship, Bentley, Rooney and Alley1986; Christianson and others, Reference Christianson2014) and mountain glaciers (e.g., Humphrey and others, Reference Humphrey, Kamb, Fahnestock and Engelhardt1993; Iverson and others, Reference Iverson, Hanson, Hooke and Jansson1995), soft-bedded processes remain poorly constrained (Walter and others, Reference Walter, Chaput and Luthi2014; Alley and others, Reference Alley, Cuffey and Zoet2019). Water flow in a soft-bedded subglacial environment has been hypothesised to occur via: Darcian flow through permeable sediments (Clarke, Reference Clarke1987); sheet flow at the ice-sediment interface (e.g., Weertman, Reference Weertman1970; Alley and others, Reference Alley, Blankenship, Rooney and Bentley1986; Flowers and Clarke, Reference Flowers and Clarke2002; Creyts and Schoof, Reference Creyts and Schoof2009); and concentrated flow in channels cut into the ice, and canals eroded into the sediment (Walder and Fowler, Reference Walder and Fowler1994; Ng, Reference Ng2000). Drainage through gaps opened and closed dynamically at the ice-sediment interface by turbulent water flow at high pressure has also been proposed as an explanation for the rapid drainage of boreholes (Engelhardt and Kamb, Reference Engelhardt and Kamb1997; Kamb, Reference Kamb2001) and both supra- and pro-glacial lakes (Sugiyama and others, Reference Sugiyama, Bauder, Huss, Riesen and Funk2008; Tsai and Rice, Reference Tsai and Rice2010, Reference Tsai and Rice2012; Hewitt and others, Reference Hewitt, Chini and Neufeld2018). Direct evidence for gap-opening at the ice-sediment interface is limited to three observational studies (Engelhardt and Kamb, Reference Engelhardt and Kamb1997; Lüthi, Reference Lüthi1999; Iverson and others, Reference Iverson2007). However, despite support from detailed analytical modelling (Schoof and others, Reference Schoof, Hewitt and Werder2012; Rada and Schoof, Reference Rada and Schoof2018) dynamic gap opening has yet to be fully developed for larger-scale numerical models of subglacial hydrology.
The water-saturated sediment layer beneath a soft-bedded ice mass can be approximated as an aquifer confined by an overlying ice aquiclude (e.g., Lingle and Brown, Reference Lingle and Brown1987; Stone and Clarke, Reference Stone and Clarke1993). And, with careful adaptation, standard hydrogeological techniques can be used to estimate subglacial aquifer properties such as transmissivity, conductivity, diffusivity and storativity. These include slug tests, where the borehole water level is perturbed by the insertion and sudden removal of a sealed pipe of known volume (Hodge, Reference Hodge1979; Stone and Clarke, Reference Stone and Clarke1993; Iken and others, Reference Iken, Fabri and Funk1996; Kulessa and Hubbard, Reference Kulessa and Hubbard1997; Stone and others, Reference Stone, Clarke and Ellis1997; Kulessa and Murray, Reference Kulessa and Murray2003; Kulessa and others, Reference Kulessa, Hubbard, Williamson and Brown2005), packer tests where the borehole is sealed near the surface and subsequently rapidly pressurised with air (Stone and Clarke, Reference Stone and Clarke1993; Stone and others, Reference Stone, Clarke and Ellis1997), and pumping tests where the borehole hydraulic head is monitored in response to water injection or extraction (e.g., Engelhardt, Reference Engelhardt1978; Engelhardt and Kamb, Reference Engelhardt and Kamb1997; Iken and Bindschadler, Reference Iken and Bindschadler1986; Lüthi, Reference Lüthi1999). Borehole drainage on connection with the bed (hereafter ‘breakthrough’), and the recovery to equilibrium water levels have also been used to determine subglacial aquifer properties (e.g., Stone and Clarke, Reference Stone and Clarke1993; Engelhardt and Kamb, Reference Engelhardt and Kamb1997; Stone and others, Reference Stone, Clarke and Ellis1997; Lüthi, Reference Lüthi1999). During breakthrough events the water level in the initially water-full borehole either: (i) drops rapidly to a new equilibrium level some tens of metres below the surface, (ii) does not drop at all, or (iii) drops slowly, or rapidly, to a new equilibrium level after a delay of minutes to days, with the variability in response usually explained in terms of the connectivity of the subglacial drainage system (e.g., Smart, Reference Smart1996; Gordon and others, Reference Gordon2001). The hydraulic conductivity of a subglacial sediment layer has also been derived from the propagation and attenuation of diurnal subglacial water pressure waves (e.g., Hubbard and others, Reference Hubbard, Sharp, Willis, Nielsen and Smart1995), and from numerical modelling of the pressure peaks induced when pressure sensors freeze in (Waddington and Clarke, Reference Waddington and Clarke1995). To date, the application of borehole response tests to marine-terminating glaciers in Greenland is limited to a single study (Lüthi, Reference Lüthi1999), presumably due to the challenges of adapting groundwater techniques to the ice-sheet setting.
The application of hydrogeological techniques requires a number of simplifying assumptions. Many techniques are fundamentally based on Darcian flow and inherently assume that the aquifer is isotropic and homogeneous; conditions that may rarely be met in the subglacial environment. Water flow in groundwater investigations is typically slow and assumed to be Darcian. While this may hold for low-velocity water flow through subglacial sediments, the discharge rates during borehole breakthrough events mean turbulent flow is likely in the vicinity of the borehole base (e.g., Stone and Clarke, Reference Stone and Clarke1993). Further complications arise due to the greater density of water than ice, overpressurising the ice at the base of water-filled glacier boreholes with the potential to raise the ice from its substrate permitting water to flow through the gap created. (Overpressure here being water pressure in excess of the ice overburden pressure). Previous studies have attempted to determine the widths of such gaps (Weertman, Reference Weertman1970; Engelhardt and Kamb, Reference Engelhardt and Kamb1997; Lüthi, Reference Lüthi1999; Iverson and others, Reference Iverson2007).
Ice boreholes provide direct access to the subglacial environment allowing sensor installation and borehole response tests. Here, we analyse borehole response tests conducted on Sermeq Kujalleq (Store Glacier) in West Greenland during summer 2019. The response tests included breakthrough events, which occurred consistently when boreholes intersected the ice-sediment interface, constant-rate pumping tests undertaken as water was pumped into the borehole as the drill stem was raised to the surface, and recovery tests following removal of the stem. The results provide insights into subglacial hydrological conditions and permit estimation of the hydraulic transmissivity and conductivity of the subglacial drainage system.
2. Methods
2.1 Field site
Sermeq Kujalleq (Store Glacier) is a major fast-moving outlet glacier of the Greenland Ice Sheet draining an ~ 34,000 km 2 catchment area (Rignot and others, Reference Rignot, Box, Burgess and Hanna2008) into Ikerasak Fjord – a tributary of Uummannaq Fjord. (Note that as several glaciers share the same name – and for continuity with previous literature – we give the English glacier name in brackets after the official Greenlandic name.) In summer 2019, we used pressurised hot water to drill seven boreholes on Sermeq Kujalleq (Store Glacier) at site R30 (N70° 34.0’, W050° 5.2’) located in the centre of the drained bed of supraglacial lake L028 (Fig. 1a; Table S1). R30 lies 30 km from the calving front at 863 m a.s.l. and is within the ablation area; there was no winter snow or firn present during the drilling campaign. Ice flow measured by a Global Navigation Satellite System (GNSS) receiver averaged 521 m a−1 in the SSW direction (217° True) between 9 July and 16 September 2019. The surface slope was calculated as 1.0° from linear regression of the ArcticDEM digital elevation model (Porter and others, Reference Porter2018) over a distance of ten ice thicknesses (10 km). Lake L028 drained via hydraulic fracture on 31 May 2019 (Chudley and others, Reference Chudley2019b) forming two major moulins (each of diameter ~6 m) located within 200 m of the drill site (Fig. 1b). Borehole-based Distributed Acoustic Sensing (DAS) in BH19c provides evidence for up to 37 m of consolidated subglacial sediment at R30 (Booth and others, Reference Booth2020), while seismic reflection surveys at site S30 (8 km to the south-east of R30; Fig. 1a) revealed up to 45 m of unconsolidated sediment overlying consolidated sediment (Hofstede and others, Reference Hofstede2018). Borehole-based investigations of englacial and basal conditions at S30 reported low effective pressures (180 − 280 kPa), an absent or thin (< 10 m) basal temperate ice layer, and internal deformation concentrated within the lowermost 100 m of ice, below the transition between interglacial (Holocene) and last-glacial (Wisconsin) ice (LGIT; Doyle and others, Reference Doyle2018; Young and others, Reference Young2019). At R30, Distributed Temperature Sensing (DTS) reveals a 70 m-thick basal temperate ice layer, the LGIT at 889 m depth and a steeply curving temperature profile with a minimum ice temperature of − 21.1°C near the centre of the ice column (Law and others, Reference Law2021).
Maps of the field site. (a) Location of the study site R30 on Sermeq Kujalleq (Store Glacier) with the location of the R29 and S30 drill sites also marked. The background is a Sentinel-2 image acquired on 1 June 2019 and the red square on the inset map shows the location in Greenland. (b) Close up of the R30 study site showing the location of boreholes, moulins and the GNSS receiver. Three boreholes intersected the ice-sediment interface (filled, colour-coded circles) and four terminated above the base (hollow circles). The background orthophoto was acquired by an uncrewed aerial vehicle survey following Chudley and others (Reference Chudley, Christoffersen, Doyle, Abellan and Snooke2019) on 21 July 2019.
2.2 Hot water drilling
Boreholes were drilled using a hot water drill system similar to that described in Makinson and Anker (Reference Makinson and Anker2014). Pressurised, hot water (11.0 MPa; ~ 80°C) was provided by five pressure-heater units (Kärcher HDS1000DE) at a regulated flow rate of 75 l min−1, through a 1,350 m long, 19.3 mm (0.75’) bore hose. A load cell and rotary encoder recorded the load on the drill tower and the hose length below the surface at 0.5 Hz with a resolution of 1 kg and 0.1 m, respectively (Figs S1–S3). Borehole logging to a depth of 325 m indicates that the hot water drilling system consistently drills boreholes that are within 1° of vertical (Hubbard and others, Reference Hubbard2021).
Boreholes (BH) were named by year and by letter in chronological order of drilling, with BH19a the first borehole drilled in 2019 (Table S1). Boreholes were drilled in two clusters with the first (BH19a, b, c and d) separated from the second (BH19e, f and g) by 70 m (Fig. 1b). Seven boreholes were drilled in 2019 with three reaching the ice-sediment interface at depths of 1,043 m (BH19c), 1,022 m (BH19e) and 1,039 m (BH19g), giving a mean ice thickness of 1,035 m and a mean elevation of the glacier sole of − 172 m a.s.l. (Table 1). Four boreholes were terminated above the ice-sediment interface (see Table S1). Prior to breakthrough boreholes were water-filled to the bare ice surface, with excess water supplied by the pressure-heater units overflowing from the top of the borehole.
Key data for the boreholes that reached the bed. Variables h 0, p w, and N were calculated for the reference period 36–60 h after each respective breakthrough, which was deemed representative of subglacial water pressure. A list of symbols is presented in Appendix A.
a Drill-indicated depths do not account for the elastic extension of the hose under load.
b Recorded in BH19e due to freeze-in of pressure transducer in BH19g.
To reduce overall drilling duration and produce a more uniform borehole radius (0.06 m four hours after termination of drilling), we optimised drilling speed using the numerical borehole model of Greenler and others (Reference Greenler2014). The borehole model was constrained by ice temperature from site R29, 1.1 km distant (Hubbard and others, Reference Hubbard2021, Fig. 1a), and a hose thermal conductivity of 0.24 W m−1 K−1. Borehole radius at the time of breakthrough was then estimated by re-running the model with the recorded drill speeds and the equilibrated ice temperature profile measured in BH19c at site R30 (Law and others, Reference Law2021). The mean borehole radius for BH19c, BH19e and BH19g output by the model at the time of borehole breakthrough was 0.07 m, with larger radii (mean of 0.10 m) in the lowermost 100 m of the ice column (Table A1) due to intentionally slower drilling as the drill approached the ice-sediment interface, together with the presence of temperate ice that was unaccounted for during initial model runs. The borehole model underestimated the near-surface (i.e. 0 − 100 m) borehole radius (r s), possibly due to turbulent heat exchange that is not included in the model, so we use the radius at the water line calculated for BH19g (0.14 m) as r s for all the borehole response tests (see Appendix B).
Analysis of the temperature time series recorded by DTS in BH19c (Law and others, Reference Law2021) shows that the boreholes rapidly froze shut. At 580 m depth, where the undisturbed ice temperature was − 21.1°C, the temperature fell below the pressure-dependent melting temperature 3 h after drilling. Within warmer ice, refreezing was slower: at 920 m depth in BH19c, the ice temperature was − 3°C and refreezing was complete after 5 d.
2.3 Pressure measurements
Basal water pressures were recorded by vibrating wire piezometers (Geokon 4500SH) installed at the base of BH19c and BH19e and a current loop transducer (Omega Engineering Ltd. PXM319) installed at the base of BH19g. Pressure records from the Geokon 4500SH were zeroed with atmospheric pressure at the surface, temperature compensated using a high-accuracy thermistor in contact with the piezometer body, and calibrated using the manufacturer's second-order polynomial to an accuracy of ± 3 kPa, equivalent to ± 0.3 m of hydraulic head. The pressure record from the PXM319 current loop transducer (accuracy = ± 35 kPa, equivalent to ± 3.6 m of head) was calibrated using the manufacturer's linear calibration and zeroed with atmospheric pressure at the surface. A pressure spike indicates that the ice surrounding the transducer installed in BH19g froze at 13.7 h post-breakthrough.
All pressure sensors were lowered until contact with the ice-bed interface was confirmed by the pressure ceasing to increase. The sensor was then raised slightly (piezometer offset: 0.05 − 0.4 m; Table 1) to prevent the piezometer from being dragged through the substrate. The borehole water level below the surface (that is the length of the uppermost air-filled section of the borehole) at installation was measured with a well depth meter, and by reference to distance markers on the piezometer cable. The final installation depth was determined by adding this water level to the depth recorded by the piezometer. The ice thickness (H i) was calculated by adding the piezometer offset to the final installation depth. Borehole positions were surveyed on 22 July 2019 using a Trimble R9s GNSS receiver with 8 min long observations post-processed using the precise point positioning service provided by Natural Resources Canada (CSRS-PPP). Borehole surface elevation was converted to orthometric EGM96 geoid heights. To allow inter-comparison of pressure records from sensors installed at different depths below the surface, water pressure was expressed as hydraulic head h, which represents the theoretical orthometric height of the borehole water level,
where ρ w = 999.8 kg m−3 is water density at 0°C, g = 9.81 m s−2 is gravitational acceleration and z is the orthometric height of the piezometer determined by subtracting the piezometer depth below the surface from the orthometric height of the borehole at the surface. Pressure was also expressed as the effective pressure N = p i − p w and the overpressure (p w − p i), the latter in respect of the excess pressure exerted at the base of water-filled boreholes due to the greater density of water than ice (Table 1). The ice-overburden pressure p i was approximated for an inclined, parallel-sided slab of ice as
where ρ i is the density of ice, H i is the height of the overlying ice column, α = 1.0° is the mean surface and bed slope (see Section 2.1), and ice density was taken as ρ i = 910 ± 10 kg m−3.
2.4 Temperature measurements
Temperature was measured using high-accuracy (± 0.05°C) thermistors (Littelfuse: PR502J2) at ~0, 1, 3, 5 and 10 m above the bed in BH19c and BH19e and also throughout the full ice column in BH19c using fibre-optic DTS (Law and others, Reference Law2021). Here we present temperature measurements recorded by the lowermost thermistor in BH19c, which was mounted with the Geokon 4500SH piezometer. We calculated the pressure-dependent melting temperature
where γ = 9.14 × 10−8 K Pa−1 is the Clausius–Clapeyron gradient determined from the basal temperature gradient (Law and others, Reference Law2021), and T tr = 273.16 K and p tr = 611.73 Pa are the triple point temperature and pressure of water, respectively.
2.5 GNSS measurements of ice motion
Time series of horizontal and vertical ice motion were determined from dual frequency (L1 + L2) GNSS data recorded by a Trimble R7 receiver at 0.1 Hz and post-processed kinematically using Precise Point Positioning with Ambiguity Resolution (CSRS PPP-AR). The GNSS antenna was mounted on a 5 m long pole drilled 4 m into the ice surface at a location between the two clusters of boreholes (Fig. 1b). Rapid re-freezing of the hole ensured effective coupling of the antenna pole with the ice. Small gaps (< 5 min) in the position record were interpolated linearly before a 6 h low pass Butterworth filter was applied. The filtered position record was differentiated to calculate velocity. The time series was then resampled to 10 min medians and a further 6 h moving average was applied to the velocity record. To prevent a shift in phase, phase-preserving filters and differentiation were used.
3. Borehole response tests
We analysed the response of borehole water pressure to the perturbations induced at breakthrough, during the continued pumping of water into the borehole while the drill stem and hose were raised to the surface, and also during the recovery phase after which borehole water pressure was in equilibrium with the pressure in the subglacial drainage system. These tests were conducted at different times since breakthrough, allowing us to investigate whether hydraulic transmissivity changed as water pressure returned to equilibrium. Rapid borehole refreezing precluded slug testing. Below we describe the borehole response test results alongside the methods.
For the majority of tests the monitoring borehole was the same as the injection borehole and these are referred to simply by the borehole name. To distinguish response tests where the injection and monitoring boreholes were different, we give the injection borehole in full followed by the monitoring borehole's letter code in brackets. A conceptual illustration of our borehole response tests is presented in Figure 2.
Conceptual diagram and nomenclature for borehole drainage via radial Darcian flow through a subglacial sediment aquifer confined by an overlying ice aquiclude. Note that monitoring boreholes are likely to have refrozen at the time of the tests and h is therefore the equivalent hydraulic head for the subglacial water pressure recorded.
All data loggers, including that of the drill, were synchronised precisely with Global Positioning System Time (GPST) immediately prior to drilling. Water pressure data were logged by separate Campbell Scientific CR1000X data loggers for each cluster of boreholes. The sampling frequency was increased to 0.2 Hz prior to borehole breakthrough, necessitating temporary suspension of thermistor measurements. Hence, no measurements of basal water temperature were made when drilling was taking place.
As it is difficult to measure the background hydraulic head without disturbing the subglacial environment it is necessary to define a reference head (h 0). The head in BH19e averaged from 36 − 60 h after BH19g breakthrough had recovered to within 0.1 m of the mean head over the 24 h period preceding BH19g breakthrough (Fig. 3b). On this basis, we define h 0 as the mean head from 36 − 60 h post-breakthrough for all tests. No corrections for background trends in hydraulic head were made but such trends are small relative to the perturbations induced (Fig. 3a).
(a) Time series of hydraulic head (h). Borehole breakthrough times are marked with a vertical dashed line and arrow. (b) Time series of head above the reference head (s = h − h 0) plotted against time since respective breakthrough for all breakthrough tests. The yellow shade marks the 24 h period selected to define h 0 (36 − 60 h post-breakthrough).
3.1 Breakthrough tests
3.1.1 Observations
All three boreholes drilled to the bed in 2019 drained rapidly upon intersecting the basal interface. During breakthrough, water levels dropped to an initial level measured during pressure transducer installation of 78, 73 and 80 m below the surface in BH19c, BH19e and BH19g (Table 1). The frictional drag of water flowing past the hose during breakthrough events caused transient ~2 kN magnitude peak forces, as recorded on the drill tower (Figs 4, S1–S3). Following the peak, force on the drill tower became constant at ~200 s post-breakthrough but at a higher level than recorded prior to breakthrough. The offset in the pre- and post-breakthrough force on the drill tower represents the difference between the weight of the hose in a water-filled and part-filled borehole.
(a) Force on the drill tower with best fit plotted against time since BH19g breakthrough, together with measured and modelled hydraulic head. (b) Volumetric flux into the subglacial drainage system (Q o) with error bars, and hydraulic head in BH19g determined by inverting the force on the drill tower. Labels (a–c) are described in Section 4.1.
As the drill stem was raised to the surface over ~2 h water continued to be pumped into the borehole, supplying an additional ~10 m3 of water (Table 1). The volume of water drained during the breakthrough events was determined from the initial water level and annular cross-sectional area of the borehole of near surface radius (r s) containing the hose of external radius (r d), yielding a mean volume for the three breakthrough events of 4.70 m3 (Table 1). Taking the duration of rapid drainage as the duration of the peak in force of ~200 s gives a mean discharge for the three breakthrough events of 2.3 × 10−2 m3 s−1 supplied from the borehole, with an additional flux supplied by the pumps Q i = 75 l min−1 (1.25 × 10−3 m3 s−1) bringing the total discharge to Q o = 2.5 × 10−2 m3 s−1, and the total volume over the ~200 s duration to 4.95 m3. The Reynolds number for outflow from the base of the borehole can be approximated as flow through a uniform cylindrical pipe, with a radius equal to that at the borehole base, the mean of which was r 0 = 0.10 m for the three boreholes (Table A1),
where η w = 0.0018 Pa s is the water viscosity at 0°C. Water flow through the boreholes near the base was turbulent with a high Re ≈ 87,500 greatly exceeding the threshold for laminar flow of 2,000 (de Marsily, Reference de Marsily1986).
3.1.2 Determining the BH19g breakthrough flux
To avoid sensor cables becoming tangled around the drill hose, pressure transducers were installed after the drill stem and hose had been recovered to the surface. Hence, no measurements of pressure were made within boreholes being drilled including during breakthrough. As the pressure response to BH19g breakthrough was captured by transducers already installed in BH19c and BH19e (Fig. 4) we now focus on the BH19g breakthrough.
We determined the time varying flux of water into the subglacial drainage system during the breakthrough of BH19g by inverting the recorded force on the drill tower from the hose, which is a combination of its weight, both in air and in water, and the frictional drag on the hose when the water drains through the borehole,
where r d is the radius of the drill, $\overline {\rho _{\rm d}}$
is the mean density of the drill (including the water core), $\Delta \overline {\rho } = \overline {\rho _{\rm d}}-\rho _{\rm w}$
, f D is the coefficient of frictional drag exerted on the outside of the hose by the down-rushing water in the borehole, H w(t) is the height of water in the borehole, F s is the force exerted by the weight of the drill stem in water, and the bulk velocity of water in the borehole during the drainage event is U w(t) = dH w/dt.
The force on the drill hose is initially set by the water height, which for a borehole full to the surface is equal to the ice thickness, therefore H w(t = 0) = H w0 = H i = 1039 m (Table 1). Since the initial force just before breakthrough F 0 = 893 N the density difference between the hose and water is
Taking ρ w = 999.8 kg m−3 gives a mean density of the hose filled with water $\overline {\rho _{\rm d}} = 1096$
kg m−3. Note that the composite density of the hose is
where ρ d is the density of the hose material, and rd = 9.7 mm is the internal bore radius of the hose. Using the calculated value of $\overline {\rho _{\rm d}} = 1096$
kg m−3 gives an estimate of the hose material density of ρ d = 1166 kg m−3, which is slightly larger than the nominal manufacturer's specification of 1149 kg m−3. This apparent extra density corresponds to an extra force measured on the drill tower prior to breakthrough of 65 N, which we interpret as a drag of 0.0625 N per metre of hose from the pumped water flowing down the centre of the hose.
Neglecting minor residual oscillations, the force F ∞ = F(t → ∞) on the drill tower after the initial rapid breakthrough was again approximately constant and is given by
From this we can infer that the final height of the water level $H_{{\rm w}\infty } = 954\pm 1\, \hbox {m}$
. That is, during BH19g breakthrough the water in BH19g transiently drops H w0 − H w∞ ≈ 85 m below the surface.
Following BH19g breakthrough a portion of the water in the borehole is rapidly evacuated into the subglacial environment. We know that the water level in the borehole decreases monotonically from an initial height H 0 to a final height H ∞ and so fit the transient response with a modified exponential solution of the form
where
A fourth order polynomial was found to be the lowest order of polynomial to accurately represent the data. The flux of water from the borehole into the subglacial environment (Q o) can then be given by
where Q i = 1.25 × 10−3 m3 s−1 is the input flux from the drill. The three constants in the polynomial y(t), c i where i = 1, …, 3, along with the drag coefficient f D were estimated using non-linear regression (MATLAB: fitnlm). The resulting constants, with error estimation, are given in Table S2. From this fit (R 2 = 0.996) of the force on the drill hose the height of water, and therefore hydraulic head, in BH19g can be calculated, together with the flux into the subglacial hydrological network (Fig. 4b). This reveals that the discharge peaked at 4.5 ± 0.1 × 10−2 m3 s−1 at 38 s after breakthrough.
3.1.3 Modelling the pressure response to BH19g breakthrough
Distinct pressure perturbations, here expressed as hydraulic head, occurred in BH19c and BH19e following the breakthrough of BH19g (Fig. 4a). In BH19e, located 4.1 m from BH19g, head instantaneously decreased by 0.93 m over a 20 ± 5 s period before rising rapidly and peaking at 14.0 m above its pre-breakthrough level 130 ± 5 s post-breakthrough. Synchronously with the drop in head observed in BH19e, a 0.11 m drop in head began in BH19c.
To analyse these pressure perturbations further we modelled the propagation of water at the contact between elastic ice and poroelastic sediment during BH19g breakthrough following Hewitt and others (Reference Hewitt, Chini and Neufeld2018). The Maxwell time for the basal temperate ice at site R30 is 10–25 min, and it is therefore reasonable to assume an elastic ice rheology for the short duration (< 4 min) of breakthrough events (Appendix C). This model accounts for pressure diffusion, flexure of the ice and deformation of the sediment, and was originally developed to describe the subglacial response to a rapidly draining supraglacial lake. The original model, which is based on Darcy's law, allowed for the formation of a subglacial cavity as well as seepage through the sediment or established subglacial networks. However, for simplicity, here we do not include cavity formation and instead assume a single effective hydraulic transmissivity for subglacial water transport; and that the fluid is incompressible. The model allows the poroelastic sediment layer to deform in response to fluid flow and pressure gradients, which allows the overlying ice to flex and bend slightly as reflected in the small (0.93 m) transient head decrease preceding the large (14.0 m) head increase recorded in BH19e following BH19g breakthrough (Fig. 4a). With these features included, the model shows how an injected fluid diffuses through the subglacial environment and how this drives a propagating flexural wave in the overlying ice.
The linearised form of the model reduces to an evolution equation for the subglacial water pressure, which for consistency is here expressed as hydraulic head h
Here A 1 = T M/b and A 2 = T B, in terms of transmissivity T, till stiffness (p-wave modulus) M, bending modulus of the ice B, and sediment thickness b. Here b is a fitting parameter, unconstrained by measurements of the actual sediment thickness, that represents the thickness of sediment affected by pressure diffusion. Assuming radial flow,
the associated flux of water q at radius r is
and q(r) = Q o(t) is the injection flux into the subglacial environment.
This problem can be solved numerically for any injection flux Q o(t). By entering the time-varying injection flux for BH19g breakthrough (Section 3.1.2) into Eqn (14), we predicted the response of hydraulic head at BH19e (4.1 m from the injection point of BH19g). An automated non-linear optimisation procedure (MATLAB: fitnlm) was used to determine the best-fit model parameters, yielding B = 2.75 × 109 Pa m3, M/b = 1 × 104 Pa m−1, and T = 1.46 × 10−4 m2 s−1. The prediction initially follows the data closely and it captures the initial decrease in BH19e hydraulic head as the flexural wave passes through (Fig. 4a). However, the model does not capture the subsequent development of the pressure recorded in BH19e; instead it predicts that the pressure drops off too rapidly after the first two minutes. We discuss this discrepancy further in Section 4.1.
3.2 Pumping tests
3.2.1 Observations
Following each breakthrough event, the hose was raised back to the surface over ~2 h (Table 1; Figs S1–S3), with the continued supply of water into the borehole functioning as a pumping test. We captured the pressure response at the base of BH19e to such a pumping test following the breakthrough of BH19g (Fig. 5). Although water was pumped down the hose while it was raised to the surface for all boreholes that reached the bed, no other pumping tests were useful as they occurred prior to the installation of pressure sensors. During the BH19g(e) pumping test the water pressure was measured in BH19e, 4.1 m distant (Fig. 5).
Time series of BH19e hydraulic head (red line) capturing the response to BH19g breakthrough and the injection of water as the hose was raised to the surface. Post-breakthrough the drill stem was kept stationary at the bed for 4 min 39 s (yellow shading). Linear fits during the three pumping test periods are shown with black lines. The light blue shade marks the period during which a piezometer was lowered into BH19g, and the dark blue shade marks the time the piezometer was temporarily snagged (see Section 4.1 for details). Labels (a–e) are also described in Section 4.1.
Starting 28 min after the breakthrough of BH19g the head in BH19e increased at a steady rate of 1.24 m h−1 (Fig. 5). This period of steady increase was interrupted by the temporary shutdown of the water supply when pressure-heater units were refuelled, with the linear increase in head resuming at the slightly higher rate of 1.36 m h−1. The rate of change of hydraulic head increased again to 7.40 m h−1 when the drill stem and hose rose above the borehole water level, indicating that, while the stem was below the water line, part of the water pumped into the borehole was replacing the reducing volume displaced by the hose as it was raised to the surface. We refer to these three periods of linearly increasing head as PT1, PT2 and PT3, respectively.
Discharge from the base of BH19g (Q o) was calculated by correcting the input flux Q i (1.25 × 10−3 m3 s−1) for storage within BH19g (Q s), and for the flux offsetting the decreasing water displacement caused by the hose as it was raised to the surface (Q d)
The pumping test was undertaken 9 d after the breakthrough of BH19e. Hence, we assume that storage within BH19e was negligible due to rapid borehole refreezing within cold ice that was present above a 70 m thick basal temperate layer (Law and others, 2021). We also consider storage within temperate ice to be negligible within the time span of our experiments due to its low permeability (e.g., 10−12 − 10−8 m2; Haseloff and others, Reference Haseloff, Hewitt and Katz2019). Q d was calculated as
where r d = 0.015 m is the hose radius and $\overline {U}_{\rm d}$
is the mean drill speed. For PT3, Q d = 0 as the drill stem and hose were above the borehole water level. Q s is the flux lost to storage in the injection borehole calculated from the rate of change in head dh/dt and the cross-sectional area of the borehole, which for PT1 and PT2 is annular as the hose was below the borehole water level
where r s = 0.14 m is the radius of BH19g at the surface (see Appendix B). For PT3
As the measurement of hydraulic head in BH19g did not start until after the pumping test, we assume that the rate of change of hydraulic head was the same in BH19g and BH19e.
These calculations reveal that during the pumping test the vast majority (90%) of water pumped into the borehole was discharged from the base (Table 2). Furthermore, this discharge from the borehole base (Q o) was remarkably steady, averaging 1.12 × 10−3 m3 s−1 with a standard deviation of 1.1 × 10−6 m3 s−1. It follows that the bulk velocity of the water ($\overline {U}_{\rm w} = Q_{\rm o} / \pi r_0^2$
) through the borehole near the base during all periods was also steady, averaging 3.2 × 10−2 m s−1 with a standard deviation of 3.1 × 10−5 m s−1.
Statistics for the BH19g(e) pumping test. V o is the volume of water discharged from the borehole base during the period. All other symbols are defined in the text.
* Calculated using the Thiem (Reference Thiem1906) method (Eqn (21)).
† Calculated using the analytical solution to the simplified Hewitt and others (Reference Hewitt, Chini and Neufeld2018) model (Eqn (23b)).
To test whether the outflow of borehole water during the pumping test was laminar or turbulent we calculated the Reynolds number (Re) using Eqn (4). During all periods, Re ≈ 3750, indicating that flow of water in the bottom of the borehole was turbulent during the pumping tests. If, however, we assume that water leaves the borehole through a gap of width δ the Reynolds number for flow through this gap is
where D h is the hydraulic diameter of the water film, r is the distance from the borehole, and ϕ is the areal fraction of the bed occupied by the gap (de Marsily, Reference de Marsily1986; Iken and others, Reference Iken, Fabri and Funk1996). For thin films with a large lateral extent D h can be approximated as 2δ (de Marsily, Reference de Marsily1986) and the equation can be simplified to
Using Eqn (20) and following the approach of Lüthi (Reference Lüthi1999), the transition from turbulent to laminar flow occurs at a distance of ~1 m from the borehole base for even the low value of ϕ = 0.1. Hence, water flow beyond this point can be treated as laminar allowing the application of standard hydrogeological techniques.
3.2.2 Hydraulic transmissivity according to the Thiem method
The hydraulic transmissivity (T s) of a subglacial sediment layer can be calculated by applying the Thiem (Reference Thiem1906) method to the pumping test data. The Thiem method assumes that a steady state has been reached within a vertically-confined, homogeneous, isotropic and incompressible aquifer with Darcian flow. In these limits the hydraulic transmissivity
where r = 4.1 m is the horizontal distance between the injection borehole (BH19g) and the monitoring borehole (BH19e), and s = h − h 0, is the mean hydraulic head (h) during the pumping test above the reference head (h 0). The radius of influence (R) is the distance to the theoretical point at which the hydraulic head remains unchanged at the equilibrium level (i.e. at radial distance R, h = h 0; s = 0; Fig. 2). (Note that the subscript in T s indicates that the method used assumes Darcian flow through sediment rather than through a gap at the ice-sediment interface, later denoted T g, or some combination of the two, for which we use T to represent the effective transmissivity.) The strong response of hydraulic head in BH19e to breakthrough in BH19g and the close agreement between head in these boreholes during the recovery phase (Fig. 3) indicates that the radius of influence is greater than the distance between BH19e and BH19g, which is 4.1 m at the surface. On the other hand, assuming a homogeneous, isotropic aquifer, the lack of a positive pressure peak in BH19c suggests the radius of influence is <70 m. Using the Thiem (Reference Thiem1906) Eqn 21, and reasonable R values of 10 and 70 m gives hydraulic transmissivity from (1.31 − 4.75) × 10−5 m2 s−1 (Table 2).
Although the Thiem (Reference Thiem1906) method is well established, it has limitations. The first is that the radius of influence R is difficult to interpret physically. The second is the requirement that a steady state has been reached. A third limitation in our application is that to calculate the flux of water leaving the base of the injection borehole (BH19g) we assume that the rate of change in hydraulic head is the same in BH19g as that recorded in BH19e.
3.2.3 Hydraulic transmissivity according to the Hewitt model
An alternative method to calculate the transmissivity from the pumping test data is through the application of an analytical solution to the simplified Hewitt and others (Reference Hewitt, Chini and Neufeld2018) model. During the pumping test Q o is steady, thereby permitting an asymptotic solution of Eqn (12) that, based on the monitoring borehole at radius r being sufficiently near to the injection borehole, gives
Hence, the predicted rate of change in hydraulic head at the nearby monitoring borehole is:
This expression is independent of parameters B, M and b and is sensitive only to the transmissivity. In principle this provides an alternative means of predicting T from the measured rate of change in hydraulic head during the pumping test, which avoids the limitations of the Thiem (Reference Thiem1906) method outlined in Section 3.2.2. This method (Eqn (23b)) gives estimates of T decreasing from 7.96 × 10−5 m2 s−1 during PT1, to 3.93 × 10−5 m2 s−1 during PT2, to 0.62 × 10−5 m2 s−1 during PT3 (Table 2).
3.3 Recovery tests
3.3.1 Observations
After water input to the borehole ceased, the borehole water pressure recovered to the reference head (h 0) over ~ 36 − 50 h (Fig. 3b; Table 1). The range in recovery times can be explained by the variable timing and magnitude of the diurnal cycle in subglacial water pressure (Fig. 3). The observed recovery curves were similar (Fig. 3b) suggesting spatially uniform subglacial hydrological conditions between boreholes. We analysed the early phase of the recovery by fitting an exponential decay curve (Weertman, Reference Weertman1970, Reference Weertman1972; Engelhardt and Kamb, Reference Engelhardt and Kamb1997) and the late phase using the Cooper and Jacob (Reference Cooper and Jacob1946) recovery test method. This provides two further estimates of hydraulic transmissivity: the first at 4 − 5 h post-breakthrough (early-phase), and the second at 14 − 27 h post-breakthrough (late-phase).
3.3.2 Exponential decay curve
The early phase of the recovery curve can be approximated as an exponential decay using the water-film model of Weertman (Reference Weertman1970, Reference Weertman1972):
where s 0 is the initial recharge at the time the pumps stopped, t is the time since the pumps stopped, and D is a time constant determined by log-linear fitting (Figs 6a–c). The water-film model, which is referred to as the gap-conduit model in Engelhardt and Kamb (Reference Engelhardt and Kamb1997), is based on the Hagen–Poiseuille equation and assumes laminar flow through a constant-width gap at the interface between the ice and a level, impermeable bed.
Recovery tests including (a–c) exponential fits (black) applied to the early stage of recovery curves plotted as hydraulic head above background (s) on the logarithmic y-axis against time (t); and (d–e) Cooper and Jacob (Reference Cooper and Jacob1946) recovery test linear-log fitting (black) applied to the late stage of the recovery curves plotted as residual drawdown (s ′) against the logarithm of the time ratio (t/t ′).
In the recovery curves of tests BH19c and BH19e the first part of the curve is missing due to the time taken to lower the pressure transducer to the bed after the drill stem was raised to the surface (Fig. 3a). Hence, s 0 was also treated as an unknown. In the BH19g(e) test the monitoring borehole was different from the injection borehole and the first part of the recovery curve was recorded. The initial BH19g(e) recovery curve was not, however, exponential and linear-log fitting was delayed for 5,000 s (83 min; Fig. 6c). After this delay the trend for BH19g(e) was quasi-exponential, in common with the other tests, and s 0 was again treated as an unknown for this test (Figs 6a–c). Hence, measured s 0 for BH19g(e) is 12.7 m and that calculated by fitting Eqn (24) is 10.1 m. The resulting time constant D was 18,200 s for BH19c, 25,000 s for BH19e and 23,000 s for BH19g(e). Rearranging Eqn (9) of Engelhardt and Kamb (Reference Engelhardt and Kamb1997) allows the gap width δ to be calculated from the time constant as
Furthermore, if we make the reasonable assumption of laminar flow at distances > 1 m from the borehole (Section 3.2), the transmissivity (T g) of a continuous porous medium equivalent to a gap of width δ is given by de Marsily (Reference de Marsily1986) as
Combining Eqns (25) and (26) (see Appendix D) allows T g to be calculated directly from the time constant (D)
For each test, two values of transmissivity were calculated, bracketing the radius of influence R to 10 − 70 m.
The results show that hydraulic transmissivity was an order of magnitude lower during the early recovery phase than during the pumping test, with hydraulic transmissivity spanning the range (1.8 − 3.5) × 10−6 m s−1 equivalent to gap widths of 0.16 − 0.20 mm for gaps covering the whole of the glacier bed (ϕ = 1; Table 3).
Results from the gap-conduit model (exponential fit). Gap width and the apparent hydraulic transmissivity were calculated for two values of the radius of influence (R = 10 and 70 m). Gap widths were additionally calculated for two areal fractions of the bed covered by the gap (ϕ = 0.1 and 1.0). The apparent gap transmissivity is independent of ϕ because gap cross-sectional area is a product of δ and ϕ.
3.3.3 Cooper and Jacob recovery tests
Hydraulic transmissivity can also be derived from the later stages of the recovery curve using the Cooper and Jacob (Reference Cooper and Jacob1946) recovery test method, providing information about the nature of the subglacial hydrological system as it returns to its original state. This method is based on the observation that, after a certain period of time, drawdown (or in our case drawup) within an aquifer at a given distance from a borehole decreases approximately in proportion to the logarithm of time since the discharge (or in our case recharge) began. The method assumes a non-leaky, vertically-confined aquifer of infinite lateral extent. Although the Theis (Reference Theis1935) method – on which the Cooper and Jacob (Reference Cooper and Jacob1946) method is based – requires a constant pumping rate, the method can be applied to a recovery test (i.e. after the pumps have ceased) using the principle of superposition of drawdown (e.g., de Marsily, Reference de Marsily1986; Hiscock and Bense, Reference Hiscock and Bense2014). Under this principle, pumping is assumed to continue uninterrupted while a hypothetical drawdown well is superimposed on the monitoring well from the time pumping stopped to exactly counteract the recharge from the pump. The residual drawup s′ is
where h, h 0, Q and T are as previously defined, and W(u) and W(u′) are well functions for the real and hypothetical boreholes where
and S is the storage coefficient, which cannot be determined using this method. In the previous two equations, t is time since the start of pumping, which for our tests is at breakthrough, and t′ is the time since the pumps stopped. As per the standard Cooper and Jacob (Reference Cooper and Jacob1946) method for pumping tests, for small values of u′ and large values of t′, the well functions can be approximated so that residual drawup can be estimated from the simplified equation
Hence, linear-log fitting allows hydraulic transmissivity (T s) to be calculated,
where Δs ′ is the rate of change of residual drawup with respect to the logarithmic time ratio. The Cooper and Jacob (Reference Cooper and Jacob1946) recovery test method described above has the advantage that the rate of recharge can be assumed to be constant, in contrast to that during an actual pumping test, which may vary (Hiscock and Bense, Reference Hiscock and Bense2014).
During the recovery phase, the sampling interval was increased from 5 to 300 s. Prior to application of the Cooper and Jacob (Reference Cooper and Jacob1946) recovery test method, the data were resampled to a constant 5 s interval and interpolated linearly. The data presented in Figures 6d–f extend from the time of pressure transducer installation at the bed (or in the case of BH19g the earlier time at which the pumps were stopped), to when diurnal pressure variations began. Fitting was applied to the later stages of the recovery curve where the trend in recharge versus the logarithmic time ratio was linear, as is required for this method to be appropriate. Accordingly, hydraulic transmissivity was calculated to be 3.0 × 10−6 m2 s−1, 2.2 × 10−6 m2 s−1 and 2.8 × 10−6 m2 s−1 for BH19c, BH19e and BH19g respectively.
4. Discussion
4.1 Hydraulic ice-sediment separation
The average drop in borehole water level during breakthrough indicates that the subglacial environment accommodated 4.70 m3 of water within 200 s. For all three boreholes that reached the bed, the delayed recovery to background levels over 36 − 50 h suggests that this breakthrough water and an additional ~10 m3 of water injected during the raise, could not be efficiently drained away from the immediate vicinity of the borehole's base. For example, recovery to the reference head took 45 h following the input of 13.6 m3 of water injected into BH19g at breakthrough and during the drill stem raise (Table 1; Fig. 3b) yielding a mean discharge of 8.4 × 10−5 m3 s−1. If the boreholes had intercepted a conduit with the capacity to drain the water away efficiently then the mean discharge rate would have been higher and the recovery time would have been shorter. Hence, it follows that at least some of this water must have been temporarily stored locally. We hypothesise that water was predominantly stored within a gap opened at the ice-sediment interface facilitated by the overpressure (913 ± 101 kPa; Table 1) exerted at the base of water-filled boreholes due to the greater density of water than ice. In the following analysis we constrain the geometry of this gap and investigate how the gap width changed through time.
An approximate calculation of the plausible range in gap width can be made for the BH19g breakthrough by assuming a uniform cylindrical subglacial water sheet with a radius ranging from 10 − 70 m (that is just greater than the distance to BH19e where a positive peak in pressure was observed and just less than the distance to BH19c where there was no positive peak in pressure). Under these assumptions, a gap width of 0.3 − 16.5 mm could accommodate the 5.17 m3 of water injected in 200 s after BH19g breakthrough. This range is consistent with a lack of discernible ice surface uplift in data collected by a GNSS receiver at R30, confirming that surface uplift was below the precision of the GNSS data of ± 50 mm (Fig. S4). Assuming a straight-sided cylinder with a volume equal to that injected during BH19g of 5.17 m3 the upper bound on the surface uplift of 50 mm provides a lower bound on the radius of the uplift of ~ 5.7 m.
Further estimates of gap widths can be determined from the hydraulic transmissivity measurements. If we assume laminar flow, which is reasonable at distances > 1 m from the borehole (see Section 3.2), the gap width (δ), equivalent to a continuous porous medium with an effective hydraulic transmissivity (T g), is given by rearranging Eqn (26)
Assuming the gap is uniformly distributed across the bed (ϕ = 1) these estimates show a decrease from 0.69 mm during breakthrough to a mean of 0.18 mm during the late recovery phase (Table 4; Fig. 7). A comparable trend was measured by Lüthi (Reference Lüthi1999) using similar methods on Sermeq Kujalleq (Jakobshavn Isbræ), with gap widths decreasing from 0.7 − 0.9 mm during a pumping test to 0.5 mm during the recovery phase. Similarly, pumping tests on a prism of simulated sediment installed beneath Engabreen yielded gap widths of 0.4 − 1.0 mm during pumping and 0.1 − 0.2 mm during recovery (Iverson and others, Reference Iverson2007). We interpret this decrease in hydraulic transmissivity and equivalent gap widths with time since breakthrough (Fig. 7) as evidence for progressive closure of gaps opened at the ice-sediment interface (in response to decreasing hydraulic head). Both our estimates, and those of Lüthi (Reference Lüthi1999) and Iverson and others (Reference Iverson2007), are lower than those of 1.4 − 2.0 mm estimated from boreholes drilled on Whillans Ice Stream (formerly Ice Stream B) in West Antarctica; however, this may, at least partly, be explained by the earlier timing made possible by measuring pressure within the Whillans boreholes while they were drilled (Engelhardt and Kamb, Reference Engelhardt and Kamb1997). The areal extent of the gap exerts a relatively weak control on gap width, with gap width approximately doubling for gaps occupying just one-tenth of the bed (ϕ = 0.1; Table 4; Fig. 7). Other lines of evidence that support the gap opening hypothesis are discussed below.
Hydraulic transmissivity (T) from multiple tests and methods plotted against time (t) since respective breakthrough. The equivalent gap width (δ) is shown on the right-hand axes for gaps covering a range of fractions of the bed (ϕ = 1 and ϕ = 0.1). Where appropriate, the range in the hydraulic transmissivity derived using radius of influence R = 10 − 70 m is shown by error bars.
The initial drop in hydraulic head in BH19e was punctuated by a 14 m increase after 20 ± 5 s, which we interpret to be the arrival of the water from the BH19g breakthrough event through a gap opened at the ice-sediment interface. The delayed arrival of the pressure increase demonstrates that no efficient hydraulic connection existed between BH19e and BH19g prior to the breakthrough of BH19g. The 20 ± 5 s delay between the start of the load increase on the drill tower and the start of the pressure increase in BH19e gives a mean velocity of the pressure pulse of 0.20 ± 0.04 m s−1. Similar pressure pulse propagation velocities of 0.08 − 0.18 m s−1 were observed on Whillans Ice Stream (Engelhardt and Kamb, Reference Engelhardt and Kamb1997). If a conduit existed between BH19g and BH19e prior to breakthrough, the pressure pulse would be transmitted at the speed of sound (1440 m s−1) and attenuated in amplitude by the viscosity of water at a rate proportional to the gap width (Engelhardt and Kamb, Reference Engelhardt and Kamb1997). The observed delay of 20 ± 5 s is four orders of magnitude longer than the expected delay of a sound wave through 4.1 m of water of 0.003 s, which confirms that no conduit existed between BH19g and BH19e prior to breakthrough. Instead, we infer that the delay represents the propagation velocity of the gap tip outwards from BH19g.
On the other hand, the disturbance in hydraulic head in BH19e caused by attempts to free a piezometer snagged at 394 m depth in BH19g, demonstrates that a hydraulic connection between the two boreholes was present at this time 2.4 h after breakthrough (Fig. 5). The piezometer in BH19g was freed after repeated pulling on the cable, which caused the hydraulic head to fluctuate in BH19e, with disturbance continuing as the piezometer was lowered to the bed. We infer that this inter-borehole transmission of pressure perturbations indicates an open gap at the ice-sediment interface at this time.
The performance of the simplified Hewitt and others (Reference Hewitt, Chini and Neufeld2018) model in predicting the pressure response to borehole breakthrough provides further evidence for gap opening. The simplified model makes a reasonable prediction of the initial pressure response in BH19e to BH19g breakthrough (Fig. 4). The model closely reproduces the small (0.93 m) drop in hydraulic head followed by the rapid rise within the first minute. This suggests that the small drop in BH19e head can be explained by the propagation of a flexural wave through the ice that is faster than the spread of water. Furthermore, the initial drop in pressure indicates that the sediment is deformable because such a drop cannot be reproduced by the model if the sediment is rigid (see Fig. 7b of Hewitt and others, Reference Hewitt, Chini and Neufeld2018). The model, however, predicts that the hydraulic head should reduce much more rapidly after the peak than was observed (Fig. 4a). Furthermore the analytical solution to the model (Eqn (23b)) predicts that ∂h/∂t should decrease non-linearly as 1/t, whereas the measured linear trends in hydraulic head during the pumping test suggest that ∂h/∂t was constant (Fig. 5). Both these disparities can be explained by gap opening.
The response of hydraulic head in BH19e to BH19g breakthrough and pumping (Figs 4 and 5) resembles the idealised pressure response of petroleum reservoirs to hydraulic fracture treatment (cf. Fig. 18a of Hubbert and Willis, Reference Hubbert and Willis1957). Specifically, the BH19g(e) breakthrough curve can be interpreted as a horizontal hydraulic fracture induced from a relatively smooth borehole, which is consistent with our interpretation of gap-opening at the ice-sediment interface induced by borehole breakthrough. We can therefore apply hydraulic fracture treatment theory to interpret the response to BH19g(e) breakthrough, as follows. After the initial drop in head, the arrival of water in BH19e is marked by a steep rise (labelled A in Figs 4a and 5), and the gradient of this increase indicates compression of the water and subglacial sediment prior to the initiation of gap opening beyond BH19e. As gap opening begins the energy stored within the compressed water and sediment is transferred to gap propagation outwards from BH19e resulting in more space for the water to occupy, and therefore lower pressure and a decrease in the gradient (dh/dt; label B). The peak in head after 130 s represents the transition from stable to unstable gap opening at the so-called ‘breakdown pressure’. The ensuing transient head decrease (label C) can be explained by the gap opening rate transiently exceeding the water input rate, and by the diffusion of unevenly distributed pressure within the immature gap. With continued water input, a steady state of gap opening was reached resulting in the linear trend in hydraulic head (label D). In our pump tests, the recharge from the pump exceeded the discharge through the gap and the borehole filled with water at a linear rate determined by the supply rate from the pumps and the extraction rate of the drill hose. That water input exceeded water output during the pumping test despite discharge rates being much lower than during breakthrough provides evidence for partial gap closure in response to reduced water pressure. When the pumps ceased, head briefly stayed constant before dropping rapidly and then transitioning into a logarithmic decay representing gap closure and reversion to Darcian flow. In petroleum engineering, the pressure at the onset of the rapid drop (label E) has been interpreted to approximate the fracture propagation pressure. For BH19g(e) this occurs at 9.290 MPa, which is comparable to the ice overburden pressure (Table 1), and is thus consistent with hydraulic ice-sediment separation. This interpretation suggests that the application of hydraulic fracture models to borehole breakthrough and pumping tests would be an improvement over hydrogeological techniques such as the Thiem (Reference Thiem1906) method, which inherently assume Darcian flow through an incompressible, isotropic aquifer. Such assumptions are unlikely to be valid if gap opening is taking place and this may explain the difference between the Thiem (Reference Thiem1906) and (analytical) Hewitt and others (Reference Hewitt, Chini and Neufeld2018) estimates of transmissivity during the pumping test (Table 4; Fig. 7).
The observation of an instantaneous drop in hydraulic head of 0.11 m in BH19c in response to BH19g breakthrough without a subsequent increase in head (Fig. 4a) also cannot be reproduced by the simplified Hewitt and others (Reference Hewitt, Chini and Neufeld2018) model; the model predicts a flexural wave that would be apparent at any fixed radius as a small pressure drop followed by a large pressure rise. We hypothesise that the drop in pressure in BH19c is caused by uplift at the BH19g injection site increasing the volume of a hydraulically-isolated cavity at BH19c, and that cavity expansion without an increase in water mass leads to a reduction in water density and pressure – that is a rarefaction. The simplified Hewitt and others (Reference Hewitt, Chini and Neufeld2018) model cannot reproduce rarefactions caused by stress transfer through the ice because it assumes that water compressibility is zero and, more fundamentally, it directly couples vertical displacement of the ice to the pressure in the subglacial environment, so that cavity expansion cannot occur without an increase in pressure (and vice versa). Further evidence for hydraulic isolation of the BH19c cavity is provided by diurnal water pressure variations that are anti-correlated with those in BH19e and ice velocity (Figs 8a.b; e.g., Murray and Clarke, Reference Murray and Clarke1995; Meierbachtol and others, Reference Meierbachtol, Harper, Humphrey and Wright2016; Lefeuvre and others, Reference Lefeuvre2018). The inference of BH19c cavity isolation is also supported by the observation that diurnal pressure variations in BH19c are manifested as small (~ 0.05° C peak-to-peak) temperature cycles recorded at the base of BH19c (Fig. 8). This demonstrates that the water temperature quickly equilibrates with the pressure-dependent ice temperature, which would occur within an isolated cavity but not in a connected conduit. We would expect that within a connected conduit a throughput of water from different regions of the bed at variable pressures and temperatures would mask the small pressure-driven diurnal variations in temperature.
Time series of (a) horizontal ice velocity, (b) hydraulic head in BH19c and BH19e, (c) temperature at the base of BH19c, and (d) pressure-dependent melting temperature T m calculated from the water pressure recorded in BH19c. Note that although the y-axes for (c) and (d) are offset the y-axis range is identical for both. The offset between measured temperature and T m can be explained by uncertainties in the sensor installation depths and the Clausius–Clapeyron gradient.
Rearranging the equation of state for water assuming mass is conserved and that temperature is constant, allows the pressure change to be related to the change in cavity volume
where V 0 and p w0 are the reference volume and pressure and β w = 5.1 × 10−10 Pa−1 is the compressibility of water. We can constrain the initial cavity geometry in two situations. First, the observation of no prior hydraulic connection between BH19e and BH19g, which were separated at the surface by 4.1 m, indicates the BH19e cavity was smaller than this distance. Second, the volume of water drained during BH19c breakthrough and the hose raise of 15.6 m3 provides an approximate maximum constraint on the BH19c cavity volume. These constraints are consistent with measurements of dye dilution in boreholes drilled on Isunnguata Sermia, which indicated cavity volumes of 7.6 ± 6.7 m3 (Meierbachtol and others, Reference Meierbachtol, Harper, Humphrey and Wright2016). Assuming the initial BH19c cavity volume was within the reasonable range of 0.5 − 15 m3 the small 0.11 m decrease in hydraulic head measured in BH19c located ~70 m distant can be explained by the contraction of the BH19c cavity of 0.3 − 8.2 × 10−6 m3. This demonstrates that, due to the low compressibility of water, the 0.11 m head decrease can be explained by a small cavity contraction of $5.5 \times 10^{-5} \% $
. Hence, we hypothesise that hydraulic ice-sediment separation caused by the overpressure at the base of BH19g caused elastic uplift of the BH19c cavity roof. The 0.11 m drop in BH19c head in response to BH19g breakthrough therefore provides direct evidence for the hypothesis of Murray and Clarke (Reference Murray and Clarke1995) that pressure variations in hydraulically-isolated cavities occur due to elastic displacement of the ice roof driven by perturbations in hydraulically-connected regions of the bed. We discuss this further in Section 4.3.
4.2 Hydraulic conductivity of subglacial sediments
We interpret the decrease in hydraulic transmissivity with time since breakthrough (Table 4; Fig. 7) as evidence for the closure of a gap at the ice-sediment interface that was opened by the overpressure at borehole breakthrough. It is notable that hydraulic transmissivity estimates derived using the Cooper and Jacob (Reference Cooper and Jacob1946) recovery tests were relatively constant (that is within 8 × 10−7 m2 s−1), despite the tests occurring over a wide range in time since breakthrough (14.1 − 27.2 h; Table 4; Fig. 7). Hence, these tests may be representative of Darcian flow through the sediment layer after gap closure. This suggestion is supported by the observation that the drawdown decreased logarithmically through time (Figs 6d–f) as is expected under Darcian flow, which is unlikely to be the case if gap closure was incomplete. Darcian flow through subglacial sediments was also inferred at site S30 from the initially logarithmic recovery in subglacial water electrical conductivity (EC) observed over 12 h following the dilution effect caused by drilling with low EC surface waters (Doyle and others, Reference Doyle2018).
When there is no flow through a gap at the ice-sediment interface, hydraulic transmissivity (T) is the hydraulic conductivity (K) integrated over the sediment thickness b
The sediment thickness at the borehole location has been estimated at $20_{-2}^{ + 17}$
m by fibre-optic distributed acoustic seismics in BH19c (Booth and others, Reference Booth2020). The full sediment thickness represents an upper limit for the calculation of hydraulic conductivity due to an increase in sediment compaction with depth, and the pressure-dependent depth limit to the diffusion of water from the ice-sediment interface (Tulaczyk and others, Reference Tulaczyk, Kamb and Engelhardt2000). For the range of hydraulic transmissivity from the Cooper and Jacob (Reference Cooper and Jacob1946) recovery tests of (2.2 − 3.0) × 10−6 m2 s−1 (Table 4), and a range of reasonable ‘hydraulically-active’ sediment thicknesses of 2 − 20 m, the hydraulic conductivity is (0.1 − 1.5) × 10−6 m s−1. This estimate is reasonable and within the range of hydraulic conductivities of glacial tills found in a range of settings by previous studies (Table 5). The Cooper and Jacob (Reference Cooper and Jacob1946) recovery test for BH19c was performed several hours earlier with respect to the time of breakthrough than those in BH19e and BH19g (Fig. 7) due to the earlier establishment of diurnal pressure variations in BH19c (Fig. 3b). If gap closure was still taking place, this earlier timing could explain the slightly higher transmissivity derived for BH19c. We also cannot exclude the possibility that water flow during breakthrough and pumping tests – or from previous natural subglacial water flow – winnowed fine particles from the upper layer of sediment, increasing the hydraulic conductivity of this layer (Iverson and others, Reference Iverson2007; Fischer and others, Reference Fischer, Iverson, Hanson, LeB Hooke and Jansson1998). As we cannot exclude winnowing, or be certain that the gap was fully closed, we interpret our estimates to represent an upper bound on the hydraulic conductivity of the sediment beneath this site.
Selected hydraulic conductivities of glacial sediments from the literature in ascending order. Sediments at the lower end of the scale (K ≤ 10−4 m s−1) were typically interpreted as unconsolidated sands and gravels, often associated with water flow from subglacial channels winnowing fine particles (Fischer and others, 1998).
Our inferred sediment hydraulic conductivity is two orders of magnitude higher than that determined from laboratory analysis of sediment retrieved from beneath Whillans Ice Stream (Engelhardt and others, Reference Engelhardt, Humphrey, Kamb and Fahnestock1990) and Trapridge Glacier in Canada (Murray and Clarke, Reference Murray and Clarke1995), see Table 5. A hydraulic conductivity of 10−7 − 10−6 m s−1 is, however, broadly consistent with the type of glacigenic sediment within core samples taken from Uummannaq Fjord. These core samples comprise glacimarine sediments deposited during the last glacial maxima including matrix supported diamict with angular to sub-angular clasts of basalt and granitic gneiss dispersed throughout a sandy mud matrix (Ó’Cofaigh and others, Reference Ó’Cofaigh2013).
Laboratory measurements of the hydraulic conductivity of glacial sediments, which inherently measure only Darcian flow, are typically a few orders of magnitude lower than field measurements (Table 5; Hubbard and Maltman, Reference Hubbard and Maltman2000), a disparity that could, at least partly, be explained by residual gap opening at the ice-sediment interface during borehole response tests (e.g., Fountain, Reference Fountain1994; Stone and others, Reference Stone, Clarke and Ellis1997). While in situ measurement of hydraulic conductivity of subglacial sediments appears to overestimate hydraulic conductivity under strict Darcian flow conditions, laboratory measurements provide little insight into the complexity of subglacial hydrological processes such as ice-sediment separation. Furthermore, as glacial sediment is by nature poorly sorted, with grain sizes ranging from boulders to clays, analysing samples that are large enough to be representative in laboratory experiments conducted at the scale necessary is more difficult than conducting in situ measurements (Clarke, Reference Clarke1987; Hubbard and Maltman, Reference Hubbard and Maltman2000). True subglacial water flow at this site may neither occur as entirely Darcian (laminar) flow through subglacial sediment nor exclusively through a gap at the ice-sediment interface, but rather a combination of the two. In any case, our in situ measurements represent a constraint on the effective hydraulic transmissivity that is independent of the process of water flow.
4.3 Implications for subglacial hydrology and basal motion
Subglacial water flow at glaciers underlain by porous sediment will naturally occur as laminar Darcian flow through interconnected pore spaces, although only insofar as the hydraulic transmissivity of the sediment is sufficient to accommodate the input of meltwater. With sustained inputs of water to the bed of many glaciers, from surface melt, for example, it may also be natural for a portion of that input to be stored temporarily in gaps opened at the ice-sediment interface, when water is delivered faster than it can permeate into the sediment below. The evidence presented herein demonstrates that the overpressure of a water-filled borehole can open a gap at the ice-sediment interface and need not directly intersect an active subglacial drainage system in order to drain. The delayed arrival of the pressure pulse in BH19e rules out the existence of sheet flow (Weertman, Reference Weertman1970; Alley and others, Reference Alley, Blankenship, Rooney and Bentley1986; Creyts and Schoof, Reference Creyts and Schoof2009), efficient conduits such as R-channels or canals (e.g., Röthlisberger, Reference Röthlisberger1972; Walder and Fowler, Reference Walder and Fowler1994; Ng, Reference Ng2000), and linked cavities (e.g., Kamb, Reference Kamb1987) prior to BH19g breakthrough, but supports the gap-opening theory of Engelhardt and Kamb (Reference Engelhardt and Kamb1997). We infer that prior to the breakthrough of BH19g, subglacial drainage at this location consisted exclusively of Darcian flow through subglacial sediments with a hydraulic conductivity K ≤ 10−6 m s−1.
Borehole drainage at the ice-sediment interface may be physically similar, but of lower magnitude, to that which occurs during the subglacial drainage of proglacial (Sugiyama and others, Reference Sugiyama, Bauder, Huss, Riesen and Funk2008), subglacial (e.g., Jóhannesson, Reference Jóhannesson2002) and supraglacial lakes (Tsai and Rice, Reference Tsai and Rice2010, Reference Tsai and Rice2012; Doyle and others, Reference Doyle2013; Dow and others, Reference Dow2015; Stevens and others, Reference Stevens2015; Hewitt and others, Reference Hewitt, Chini and Neufeld2018) via a broad, turbulent and transient sheet. We note that gap opening at the ice-sediment or ice-bedrock interface is conceptually the same as horizontal hydraulic fracture at this interface as envisaged by previous studies (Tsai and Rice, Reference Tsai and Rice2010, Reference Tsai and Rice2012; Hewitt and others, Reference Hewitt, Chini and Neufeld2018). Rapid water flow into this narrow gap is likely to be turbulent (Section 3.1.1); however, flow must become laminar near the gap tip as the width of the gap decreases to zero, and flow velocity will also decrease with distance from the injection point (Hewitt and others, Reference Hewitt, Chini and Neufeld2018). Continued sheet flow through a uniform gap would be unstable as irregularities in flow would theoretically favour the formation of conduits through preferential sediment erosion and concentrated ice melt from frictional heat (Röthlisberger, Reference Röthlisberger1972; Walder and Fowler, Reference Walder and Fowler1994; Ng, Reference Ng2000). Conduit development beneath kilometre-thick ice is, however, anticipated to require continuous water supply at high pressure over prolonged periods, which may only occur if there is continued water input from the surface (e.g., Dow and others, Reference Dow, Kulessa, Rutt, Doyle and Hubbard2014, Reference Dow2015). Hence, our inference of complete, or at least partial, gap closure in response to declining pressure is consistent with existing theory as the water volumes provided by borehole drainage and subsequent pumping (~15 m3) are likely insufficient to establish an efficient conduit beneath kilometre-thick ice. The development of efficient conduits in response to borehole breakthrough can also be excluded by the low discharge rate of 8.4 × 10−5 m3 s−1 calculated from the 45 h required for hydraulic head to recover to the equilibrium level following the injection of 13.6 m3 of water at BH19g breakthrough and during the drill stem raise. Although we cannot rule out the persistence of stable sheet flow following borehole drainage facilitated by clasts partially supporting the ice overburden pressure (Creyts and Schoof, Reference Creyts and Schoof2009), our observations of a progressive decrease in hydraulic transmissivity can be entirely explained by gap closure and a reversion to Darcian flow through the sediment layer. For simplicity, this and previous studies (Tsai and Rice, Reference Tsai and Rice2010, Reference Tsai and Rice2012; Hewitt and others, Reference Hewitt, Chini and Neufeld2018) make the reasonable assumption that initial gap opening is elastic; however, where temperate ice is present, as it is at R30, viscous deformation cannot be neglected during the longer time scales of pumping tests or lake drainage events (Appendix C). The application of a viscoelastic model (e.g., Reeh and others, Reference Reeh, Christensen, Mayer and Olesen2003) to borehole response tests (and lake drainage events) would therefore represent an improvement over the analysis presented herein (and previously).
The instantaneous 0.11 m drop in BH19c head in response to BH19g breakthrough (Fig. 4a) provides direct evidence for the hypothesis of Murray and Clarke (Reference Murray and Clarke1995) that pressure variations can be transmitted to unconnected cavities through elastic displacement of the ice roof. Murray and Clarke (Reference Murray and Clarke1995) theorised that uplift caused by high water pressure relieves the pressure in adjacent hydraulically-isolated cavities. This hypothesis is one of three hypotheses of mechanical forcing of water pressure that have been proposed to explain the often observed diurnal variation of water pressure in hydrologically-isolated cavities that is out of phase with both ice velocity and water pressure in boreholes and moulins deemed to be connected to efficient subglacial conduits (Murray and Clarke, Reference Murray and Clarke1995; Engelhardt and Kamb, Reference Engelhardt and Kamb1997; Gordon and others, Reference Gordon1998; Dow and others, Reference Dow, Kavanaugh, Sanders, Cuffey and MacGregor2011; Andrews and others, Reference Andrews2014; Ryser and others, Reference Ryser2014; Lefeuvre and others, Reference Lefeuvre, Jackson, Lappegard and Hagen2015; Meierbachtol and others, Reference Meierbachtol, Harper, Humphrey and Wright2016; Rada and Schoof, Reference Rada and Schoof2018). While we cannot rule out the possibility that such anti-correlated diurnal pressure and velocity variations in BH19c (Fig. 8) can be attributed to the alternative hypotheses of cavity expansion and contraction caused by longitudinal strain (Ryser and others, Reference Ryser2014) or basal sliding (Iken and Truffer, Reference Iken and Truffer1997; Bartholomaus and others, Reference Bartholomaus, Anderson and Anderson2011; Hoffman and Price, Reference Hoffman and Price2014), displacement of the ice roof due to elastic uplift during gap-opening at BH19g breakthrough can entirely explain the 0.11 m instantaneous drop in BH19c head. It is therefore plausible that elastic displacement of the ice roof by diurnal pressure variations within a nearby conduit also explains the anti-correlated diurnal variations in BH19c pressure. This assertion is supported by three-dimensional full-Stokes modelling (Lefeuvre and others, Reference Lefeuvre2018) that reproduced anti-correlated pressure variations between connected and unconnected components of the subglacial drainage system without invoking cavity expansion caused by sliding.
Similar to borehole breakthrough events, we argue that water flow at the ice-sediment interface may also occur at times of naturally high subglacial water pressures. It is important to note that the gap widths we report are probably larger than would have occurred naturally for the same volume of cold glacial water because warm drilling water would have enlarged the gaps through ice melt. The greater variability in meltwater supply means that gap opening at the ice-sediment interface is more likely to occur naturally on the Greenland Ice Sheet, and on mountain glaciers, than on the West Antarctic ice streams where the process was originally inferred (Engelhardt and Kamb, Reference Engelhardt and Kamb1997). Hence, gap opening at the ice-sediment interface has important implications for our understanding of subglacial hydrological systems that extends beyond its ability to explain the drainage of boreholes. Subglacial hydrology in ice-sheet models may for instance include exchanges of water flowing partly at the interface and partly within subglacial sediment, which has proven efficient in reproducing day-to-day variations in ice flow as observed at the land-terminating southwest ice margin (Bougamont and others, Reference Bougamont2014). Darcian flow and gap-opening therefore provide a physical explanation for the partitioning of water flowing at the interface and within subglacial sediment.
Gap-opening may also play a role in the formation and growth of subglacial drainage systems. Within the framework of existing theory, gap opening provides the initial conduit that may later develop into an inefficient narrow orifice in a distributed (i.e. linked cavity) drainage system (Kamb, Reference Kamb1987), which may ultimately develop into an efficient channel or canal (Röthlisberger, Reference Röthlisberger1972; Walder and Fowler, Reference Walder and Fowler1994; Ng, Reference Ng2000). That the overpressure of a water-filled vertical conduit stretching from the surface to the bed (i.e. a borehole) can open a gap at the ice-sediment interface, despite the low volumes of water involved, has implications for the establishment of subglacial drainage of the much larger water volumes supplied via moulins, crevasses and supraglacial lakes. It illustrates the manner in which regions of the basal environment can become hydrologically connected during peaks in water pressure. Hence, gap opening can explain transient periods of borehole water pressure synchroneity that abruptly punctuate the often observed long-term pattern of anti-correlated variations in water pressure and velocity measured in hydraulically-isolated cavities during periods of high water pressure (e.g., Murray and Clarke, Reference Murray and Clarke1995; Engelhardt and Kamb, Reference Engelhardt and Kamb1997; Harper and others, Reference Harper, Humphrey, Pfeffer and Lazar2007; Andrews and others, Reference Andrews2014; Rada and Schoof, Reference Rada and Schoof2018). If areas of the bed that were previously hydraulically isolated experience net drainage as a result of gap opening at the ice-sediment interface, it may also explain the hydro-mechanical regulation of ice flow (e.g., Sole and others, Reference Sole2013; Tedstone and others, Reference Tedstone2015; Davison and others, Reference Davison2020), which observations suggest cannot be entirely explained by water pressures within efficient channels (Andrews and others, Reference Andrews2014). It follows that drainage at the ice-sediment interface and Darcian flow through sediments with a low hydraulic conductivity may be two of potentially multiple processes behind the hypothesised weakly-connected component of the subglacial drainage system (Hoffman and others, Reference Hoffman2016).
A drainage system consisting of cavities, which we assume are present at the base of our boreholes, linked via gaps opened at the ice-sediment interface would at first appear similar to the linked cavity theory of glacial drainage, which consists of cavities connected via narrow orifices (e.g., Kamb, Reference Kamb1987). There is, however, an important distinction in that the linked cavity model specifies that orifices are continuously open and water flow is inefficient and turbulent due to the length and narrowness of orifices (Kamb, Reference Kamb1987). Modification of the linked cavity theory to allow transient gap opening between cavities under high water pressure with turbulent flow would explain the same characteristics associated with linked cavity drainage systems: enhanced basal motion, sediment entrainment (as indicated by increased turbidity) and increased connectivity of the bed at times of high water pressure. It would also explain the existence of neighbouring yet behaviourally-independent subglacial drainage subsystems in close proximity (e.g., Murray and Clarke, Reference Murray and Clarke1995; Harper and others, Reference Harper, Humphrey, Pfeffer and Lazar2007; Rada and Schoof, Reference Rada and Schoof2018), which the majority of previous models of subglacial drainage cannot reproduce as they inherently allow water to diffuse across the entire glacier bed (e.g., Schoof, Reference Schoof2010; Hewitt, Reference Hewitt2013; Werder and others, Reference Werder, Hewitt, Schoof and Flowers2013). This implies a strong link between subglacial hydrology, stresses within the ice and basal motion that will be challenging to reproduce within numerical models due to the requirement to combine linear-elastic gap opening with a viscous ice rheology.
To date, every borehole drilled on Sermeq Kujalleq (Store Glacier) drained rapidly and immediately upon reaching the bed. This includes three boreholes at R30 in 2019, four boreholes at R29 in 2018 (unpublished), and seven boreholes at S30 in 2014 and 2016 (Doyle and others, Reference Doyle2018). A similar pattern of rapid borehole drainage, with a small number of exceptions, has been reported for Whillans Ice Stream in West Antarctica (Engelhardt and Kamb, Reference Engelhardt and Kamb1997) and Sermeq Kujalleq (Jakobshavn Isbræ) in West Greenland (Lüthi, Reference Lüthi1999). While the results presented here provide further evidence for gap opening as a mechanism for rapid borehole drainage, it also raises the question of why some boreholes on other ice masses do not drain rapidly upon reaching the bed. Some boreholes appear to never drain (e.g., Smart, Reference Smart1996), while others drain slowly (e.g., Andrews and others, Reference Andrews2014), and others drain after a delay (e.g., Gordon and others, Reference Gordon2001; Kamb and Engelhardt, Reference Kamb and Engelhardt1987; Engelhardt and Kamb, Reference Engelhardt and Kamb1997; Fischer and Clarke, Reference Fischer and Clarke2001). This heterogeneity, which often occurs within the same field site, could be explained by the stress regime, boreholes terminating blind in debris-rich basal ice before they are able to connect to the subglacial drainage system, or by the presence of impermeable barriers such as areas of ice-bedrock contact or cold ice, the latter of which can occur even within predominantly temperate glaciers (Robin, Reference Robin1976). A detailed discussion of the heterogeneity of borehole drainage is not warranted here (see instead Smart, Reference Smart1996 and Gordon and others, Reference Gordon2001), but we do seek an explanation for the homogeneity in borehole drainage observed to date on Sermeq Kujalleq (Store Glacier). Hot water drilling is ineffective at penetrating debris-rich basal ice, which is characteristic of many exposed margins of the Greenland Ice Sheet, for example, on Russell Glacier (Knight and others, Reference Knight, Waller, Patterson, Jones and Robinson2002) and at the base of icebergs discharging from Sermeq Kujalleq (Jakobshavn Isbræ; Lüthi and others, Reference Lüthi, Fahnestock and Truffer2009), yet none of the boreholes drilled to date on Sermeq Kujalleq (Store Glacier) terminated above the bed due to an obstruction by englacial clasts. We therefore speculate (while noting the small number of boreholes drilled at a limited number of sites) that debris content within basal ice on Sermeq Kujalleq (Store Glacier) may be low. If so, this could be explained by the removal of debris-rich basal ice formed upstream by basal melt. Furthermore, low (and potentially even negative) effective pressures (e.g. − 46 ± 102 kPa at R30; Table 1) are conducive to hydraulic ice-bed separation (e.g., Schoof and others, Reference Schoof, Hewitt and Werder2012) and these conditions are found at all the Sermeq Kujalleq (Store Glacier) sites drilled to date. Modelling of subglacial drainage through a poroelastic sediment and cavity suggests that elastic gap opening is enabled by the suction of water from an underlying porous sediment layer without the requirement for a pre-wetted water film (Hewitt and others, Reference Hewitt, Chini and Neufeld2018). We therefore conclude that rapid borehole drainage on Sermeq Kujalleq (Store Glacier) is facilitated by low effective pressures, subglacial sediment, and a potentially low debris content within basal ice.
Booth and others (Reference Booth2020) used the low basal reflectivity in vertical seismic profiles to infer that the subglacial sediment layer at site R30 has an acoustic impedance similar to that of basal ice, and from this, they suggested that the sediment is consolidated, and neither deforming nor lithified. The inference that the sediment layer is not deforming implies that the fast ice velocity at this site must be accommodated by either enhanced internal deformation of the ice, ice-sediment decoupling under high water pressure (e.g., Iverson and others, Reference Iverson, Hanson, Hooke and Jansson1995, Reference Iverson2007), or deformation of a sediment layer thinner than the 5 − 10 m vertical resolution of the seismic technique. With regard to the last assertion we note that sediment deformation often occurs within an upper layer that is typically only decimetres to a few metres thick (e.g., Clarke, Reference Clarke1987; Murray, Reference Murray1997; Humphrey and others, Reference Humphrey, Kamb, Fahnestock and Engelhardt1993; Engelhardt and Kamb, Reference Engelhardt and Kamb1998), and that the shape of the pressure pulse during BH19g breakthrough can only be reproduced using the model of Hewitt and others (Reference Hewitt, Chini and Neufeld2018) if the sediment layer is deformable. While the extent of sediment deformation beneath this site remains inconclusive the evidence presented herein supports the hypothesis of ice-sediment decoupling under periods of high water pressure. Indeed, we suggest that the theory of gap opening at the ice-sediment interface (Engelhardt and Kamb, Reference Engelhardt and Kamb1997) may involve the same physical process as ice-sediment decoupling envisaged by Iverson and others (Reference Iverson, Hanson, Hooke and Jansson1995). To explain the reverse tilt of inclinometers just below the ice-sediment interface, Iverson and others (Reference Iverson, Hanson, Hooke and Jansson1995) envisaged that sediment would be squeezed into the zone of uplift at times of high water pressure. The modulation of slip by pressurised water at the ice-sediment interface was confirmed by pump tests on a simulated prism of till on Engabreen (Iverson and others, Reference Iverson2007). Further evidence for gap opening and decoupling at the ice-sediment interface is provided by (as far as we are aware) unrepeated, direct observation of a cm-wide gap at the ice-sediment interface of Blue Glacier, USA (Engelhardt and others, Reference Engelhardt, Harrison and Kamb1978). Borehole photography revealed a ~0.1 m thick sediment layer overlying bedrock that was mechanically and visibly distinct from a 0.1 − 16.0 m thick debris-laden basal ice layer. Engelhardt and others (Reference Engelhardt, Harrison and Kamb1978) suggested that the gap was opened by the overpressure of the water-filled borehole and that basal sliding velocities were faster where gaps were present. They also inferred that interstitial pressure within the sediment must be close to or at the ice overburden pressure in order to prevent the basal ice merging with the sediment layer through regelation, an assertion supported by Rempel (Reference Rempel2008). Hence, further in situ observations are required to investigate whether ice-sediment decoupling occurs via a gap at the ice-sediment interface or through an increase in the thickness of the sediment layer as proposed by Iverson and others (Reference Iverson, Hanson, Hooke and Jansson1995), or a combination of both processes as modelled by Hewitt and others (Reference Hewitt, Chini and Neufeld2018).
5. Conclusions
Detailed measurements of pressure pulses during a borehole breakthrough event, and a decrease in hydraulic transmissivity with time since breakthrough, provide evidence for hydraulic gap opening and closure at the ice-sediment interface, with gaps opening and closing in response to water pressure. Analysis of the subsequent recovery of subglacial water pressure indicates that the hydraulic conductivity of the subglacial sediment layer is on the order of 10−7 − 10−6 m s−1, which suggests it is coarse-grained and more permeable than the fine-grained sediments beneath West Antarctic ice streams. As seismic surveys suggest that sediment at this site is not deforming, we infer that fast basal motion may be accommodated by ice-sediment decoupling and potentially shallow-depth sediment deformation in a layer thinner than the 5 − 10 m resolution of the seismic technique.
Observations of a pressure drop simultaneous with the breakthrough of a borehole 70 m away provides direct evidence for the hypothesis that anti-correlations between water pressure in connected and unconnected regions of the bed can be explained via elastic displacement of the ice roof.
We argue that water flow via gaps opened at the ice-sediment interface is likely to play a critical role in both basal motion and the development of subglacial hydrology on soft-bedded ice masses, and that Darcian flow through sediments may explain the drainage and recharge of areas of the bed that are otherwise hydrologically isolated.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/jog.2021.121.
Data
The data sets presented in this paper are available for download from https://doi.org/10.6084/m9.figshare.16838020.
Acknowledgements
This research was funded by the European Research Council as part of the RESPONDER project under the European Union's Horizon 2020 research and innovation program (Grant 683043). T.R.C. and R.L. were supported by Natural Environment Research Council Doctoral Training Partnership Studentships (Grant NE/L002507/1). We thank Lee Greenler for providing the code for modelling borehole diameter; Katie Miles, Emma Docherty and Tom Chase for assistance in the construction of borehole sensor strings; and Sean Peters, Mickey MacKie, Mike Prior-Jones, Eliza Dawson and Tun Jan Young for assistance in the field. We are very grateful to Ann Andreasen and the Uummannaq Polar Institute for their kind hospitality. Comments from the Scientific Editor, Ali Graham, B. de Fleurian and two anonymous reviewers resulted in significant improvements to this paper.
Author contributions
The overall research project (RESPONDER) was led by P.C., with B.H. leading the hot water drilling and borehole instrumentation reported herein. Data collection was led by S.D., with contributions from B.H., P.C., R.L., C.S. and T.C. S.D. conducted the data analysis. R.L. adapted and ran the borehole drilling model. T.C. surveyed the borehole positions and led site mapping. D.H. and J.N. calculated the breakthrough volumetric flux and pressure response. The manuscript was written by S.D., with contributions from all co-authors.
Appendix A. List of symbols
- α
Surface and bed slope (°)
- β w
Water compressibility (5.1 × 10−10 Pa−1)
- γ
Clausius–Clapeyron constant (9.14 × 10−8 K Pa−1)
- δ
Gap width (m)
- η i
Effective ice viscosity (Pa s−1)
- η w
Water viscosity at 0°C (0.0018 Pa s)
- ρ i
Ice density (910 ± 10 kg m−3)
- ρ w
Water density at 0°C (999.8 kg m−3)
- ρ d
Hose density (kg m−3)
- τ e
Effective stress (Pa)
- ϕ
Areal fraction of the bed covered by gap
- A
Rate factor in Glen's flow law (Pa−3 s−1)
- b
Sediment thickness (m)
- B
Bending modulus of the ice (Pa m3)
- D
Time constant (s)
- E
Elastic modulus of ice (9.3 GPa)
- f
Shape factor
- f D
Frictional drag coefficient
- F
Force on the drill tower (N)
- g
Gravitational acceleration (9.81 m s−2)
- h
Hydraulic head (m)
- h 0
Reference hydraulic head (m)
- H i
Ice thickness (m)
- H w
Water height (m)
- K
Hydraulic conductivity (m s−1)
- M
Sediment stiffness (p-wave modulus) (Pa)
- n
Exponent in Glen's flow law (3)
- N
Effective pressure (Pa)
- p i
Ice overburden pressure (Pa)
- p w
Subglacial water pressure (Pa)
- p tr
Triple point pressure of water (611.73 Pa)
- Q
Volumetric flux (m3 s−1)
- r
Radial distance (m)
- r d
External hose radius (0.015 m)
- r 0
Borehole radius at base (m)
- r s
Borehole radius at near-surface (m)
- R
Radius of influence (m)
- Re
Reynolds number
- s
Recharge (s = h − h 0) (m)
- s 0
Reference recharge (m)
- S
Storage coefficient (m)
- t
Time (s)
- t M
Maxwell time (s)
- T
Hydraulic transmissivity (m2 s−1)
- T m
Melting temperature of ice (°C)
- T tr
Triple point temperature of water (273.16 K)
- U d
Drill velocity (m min−1)
- U w
Water velocity (m s−1)
- V
Volume (m3)
- W(u)
Well function
- z
Orthometric height (m)
Appendix B. Borehole radius
As the hose radius (r d) and speed (U d) are known, the differential rate of change in hydraulic head below and above the water line during the BH19g(e) pumping test allows the borehole radius at the water line (r s) to be determined as follows. The total volumetric flux of water stored within the borehole when the drill hose was below the water line during PT2 is Q b2 = Q s2 + Q d2, or alternatively
where the numeric subscript indicates the period. Similarly the borehole storage flux with the drill stem above the water line during PT3 is
Assuming water input (Q i) and output (Q o) were constant at the transition from PT2 to PT3
Therefore equating fluxes gives
Expanding on the left-hand side gives
Rearranging gives
and factorising gives
which we rearrange to find
Using Eqn (B8), the known hose radius (r d= 0.015 m), the measured mean drill speed during PT2 ($\overline {U}_{\rm d} = 8.82$
min−1), and the rate of change in hydraulic head during PT2 (dh 2/dt = 1.36 m h−1) and PT3 (dh 3/dt = 7.40 m h−1), gives a borehole radius at the water-line r s = 0.14 m. This estimate is double that of the borehole model (r s = 0.07 m; Table B1), but consistent with the borehole radius measured at the surface.
Measurements were not made of BH19g but BH19e had a radius at the surface of 0.17 m. As the pumping test period was not recorded in BH19c and BH19e we assume that their near-surface radius was the same as BH19g: that is, we assume r s = 0.14 m for all response tests. Near-surface borehole radii larger than predicted by the Greenler and others (Reference Greenler2014) model could be explained by turbulent heat exchange from warm upwelling water. Laminar flow is specified in the model. The effect of turbulent heat exchange on borehole radius would decrease with depth so the model should perform better near the base. With no better estimate available, we therefore use the model output for the borehole radius at the base (r 0; Table B1).
Borehole radii at the time of borehole breakthrough predicted using the model of Greenler et al. (2014) over ten depth intervals ranging from the ice surface to the ice-sediment interface at a depth below the ice surface corresponding to the ice thickness (H i)
Appendix C. Elastic response of ice to borehole breakthrough
Here we consider the relative importance of viscous and elastic deformation in the response of the ice sheet at site R30 to borehole breakthrough forcing by calculating the Maxwell relaxation time
where E = 9.3 GPa is the elastic modulus for ice (Sinha, Reference Sinha1978), and η i is the effective ice viscosity. The effective viscosity can be given as
where A and n = 3 are the rate factor and exponent in Glen's flow law, and τe is the effective stress (Hutter, Reference Hutter1983). For simplicity, we estimate the effective stress as
where, for site R30, f ≈ 0.75 is the shape factor representing the proportion of driving stress supported by basal drag (Nye, Reference Nye1952). Using Eqn (C3), the effective stress at site R30 is 121 kPa. We assume that viscous deformation will be greatest within the basal temperate ice layer and therefore use upper and lower limits of A for temperate ice of 5.5 − 2.4 × 10−24 Pa−3 s−1 (Cuffey and Paterson, Reference Cuffey and Paterson2010). With these values the effective viscosity is 6.2 − 14.2 × 10−12 Pa s−1, and the Maxwell time is 11 − 25 min. Hence, assuming elastic ice rheology at site R30 is reasonable during the initial stages of gap opening. Over the time scales relevant to pumping and recovery tests viscous deformation should not be neglected and a viscoelastic model (e.g., Reeh and others, Reference Reeh, Christensen, Mayer and Olesen2003) would be more appropriate. Note, however, that the rheology of the ice actually drops out of the asymptotic solution of the Hewitt and others (Reference Hewitt, Chini and Neufeld2018) model in Eqn (22), and so incorporating viscous deformation may not have a large effect on the predictions of transmissivity from that model.
Appendix D. Transmissivity from time constant
The hydraulic transmissivity (T g) of a porous medium equivalent to a gap of uniform width δ is given by de Marsily (Reference de Marsily1986) as
The time constant D is given by
which is Eqn (7a) of Weertman (Reference Weertman1970) and Eqn (9) of Engelhardt and Kamb (Reference Engelhardt and Kamb1997). Combining Eqns (D1) and (D2) as follows allows the hydraulic transmissivity to be approximated from the time constant D. Inserting ϕ and then multiplying both sides of Eqn (D2) by two gives
This permits simplification by inserting the inverse of Eqn (D1) into Eqn (D3)
Multiplying both sides by T g gives
And further rearranging gives
which is Eqn (8.7) of Lüthi (Reference Lüthi1999) and Eqn (27) of this paper.
Open access
