Nepal, revealed by dye tracing Surface and subsurface hydrology of debris-covered Khumbu Glacier, Earth and Planetary Science Letters

Contents lists available atScienceDirect

Earth and Planetary Science Letters

www.elsevier.com/locate/epsl

Surface and subsurface hydrology of debris-covered Khumbu Glacier, Nepal, revealed by dye tracing

Katie E. Miles

, Bryn Hubbard

, Duncan J. Quincey

, Evan S. Miles

, Tristram D.L. Irvine-Fynn

, Ann V. Rowan

aCentreforGlaciology,DepartmentofGeographyandEarthSciences,AberystwythUniversity,Aberystwyth,SY233DB,UK bSchoolofGeography,UniversityofLeeds,Leeds,LS29JT,UK

cSwissFederalResearchInstituteWSL,8903Birmensdorf,Switzerland dDepartmentofGeography,UniversityofSheffield,Sheffield,S102TN,UK

a r t i c l e i n f o a b s t ra c t

Articlehistory:

Received29October2018

Receivedinrevisedform14February2019 Accepted17February2019

Availableonline6March2019 Editor: A.Yin

Keywords:

glacierhydrology dyetracing debris-coveredglacier Himalaya

subglacial

While the supraglacial hydrology ofdebris-covered glaciers is relatively well studied,almost nothing is known about how water is transported beneath the glacier surface. Here, we report the results of sixteen fluorescent dye tracing experiments conducted in April–May 2018 over the lowermost 7 km ofthehigh-elevation,debris-coveredKhumbuGlacier,Nepal,tocharacterisetheglacier’s surface and subsurface drainagesystem. Dyebreakthroughs indicatedalikely highlysinuous and channelised subsurface hydrologicalsystemdrainingwaterfromtheupper partoftheablation area.Thisflowpath was distinct from the linked chain of supraglacial ponds present along much of the glacier’s lower ablationarea,throughwhichwaterflowwasextremelyslow(∼^{0.003m s−1),likelyreflectingthestudy’s} timing duringthe pre-monsoon period.Subsurface drainagepathways emerged at theglacier surface close to the terminus, and flowed intosmall near-surface englacial reservoirs that typically delayed meltwatertransitbyseveralhours.Weobservedrapidpathwaychangesresultingfromsurfacecollapse, indicatingafurtherdistinctiveaspectofthedrainageofdebris-coveredglaciers. Weconcludethatthe surface and subsurface drainageof Khumbu Glacieris bothdistinctive and dynamic, and argue that furtherinvestigationisneededtorefinethecharacterisationandtestitsregionalapplicabilitytobetter understandfutureHimalayandebris-coveredglaciermeltwaterdeliverytodownstreamareas.

©^{2019TheAuthor(s).PublishedbyElsevierB.V.ThisisanopenaccessarticleundertheCCBYlicense} (http://creativecommons.org/licenses/by/4.0/).

1. Introduction

Meltwater from Himalayan glaciers and snow feeds some of Earth’s largest river systems, influencing the supply of water to

∼^{1.4 billionpeople (Barnett etal.,} 2005; Bolch, 2017; Immerzeel et al., 2010). Approximately 30% of Himalayan glaciers have a supraglacialdebriscover(Thakurietal.,2014) thatinfluencesmass lossprocesses,thelociofmaximummelt,andthevolumeofmelt- waterproduced (Luckman etal., 2007; Østrem, 1959; Thompson etal., 2016). Where thedebris coverexceeds severalcentimetres in depth andcovers a considerable portion of the glacier’s abla- tionarea,itinfluencesthehydrologicalsystembothonthesurface andbelowit(Fyffeetal.,2019).Theseextensivedebriscoverspro- ducearangeofsurfacefeaturesnotcommonlyfoundonclean-ice glaciers,suchassupraglacialpondsandicecliffs.

*

E-mailaddress:kam64@aber.ac.uk(K.E. Miles).

The supraglacial hydrology of debris-covered glaciers has re- ceivedincreasingattentioninrecentyears.Forexample,itiswell- documented that the formation of proglacial moraine-dammed lakes, and the consequent presence of a local base-level, can be facilitated by the development, growth, and coalescence of supraglacialpondsintoachainoflinkedponds(Bennetal.,2012;

Mertes et al., 2016; Sakai, 2012; Thompson et al., 2012). Melt rates are disproportionately high at pond margins due to the continued horizontal and vertical incision (Miles et al., 2016;

Sakaietal.,2000),andrecentstudieshavefoundthatsupraglacial pondsareexpandingtocoveran increasingproportionofthesur- face of debris-covered glaciers (Gardelle et al., 2012; Watson et al., 2016).This hasimplications forgreatermeltwater production and water storage, since ponds moderate diurnal glacier runoff (Irvine-Fynnetal.,2017).

On debris-covered glaciers, supraglacial streams do not tend to persistforlong distances,instead incisingtobecome englacial features (Gulley etal., 2009; Iwata etal., 1980). Supraglacialand shallow englacial conduits located towards the centre of debris- covered glaciers havebeen suggestedto bediscontinuous dueto https://doi.org/10.1016/j.epsl.2019.02.020

0012-821X/©^{2019TheAuthor(s).PublishedbyElsevierB.V.ThisisanopenaccessarticleundertheCCBYlicense}(http://creativecommons.org/licenses/by/4.0/).

thevariable surface topography, andare therefore mostlikely to transportmeltwaterbetweensupraglacialponds(Bennetal.,2017;

Miles et al., 2017a; Narama et al., 2017; Thakuri et al., 2016).

Longer-distance transport has been observed through perennial sub-marginalchannels(Bennetal.,2017; Thompsonetal.,2016), while more conduits open and transport water more effectively throughthemonsoonal meltseason thanthrough thedryseason (Bennetal.,2012; Hewittetal.,1989;Milesetal.,2017a,2017b;

Sakaietal.,2000).

Aside from a limited number of speleological investigations (Benn et al., 2017; Gulley et al., 2009; Gulley and Benn, 2007;

Narama et al., 2017), most work on Himalayan debris-covered glacier englacial drainage has been inferred from proxies. Even lessis known about possible subglacial drainage networks; such systems have been deduced on the basis of methods such as proglacial sediment and water analyses (Haritashya et al., 2010;

Hasnainand Thayyen, 1994) and remote-sensing observations of seasonal changes in glacier surface velocity (Benn et al., 2017;

Copland et al., 2009; Kraaijenbrink et al., 2016; Quincey et al., 2009). Seasonal systemevolution has, however,been determined directly from a few dye-tracing studies (Hasnain et al., 2001;

Liuetal.,2018; Pottakkaletal.,2014).Fyffeetal. (2019) alsocon- ductedadyetracingstudyonadebris-coveredglacierintheItal- ianAlps(MiageGlacier),reportingthatthecontinuousdebriscover intheglacier’slowerablationareaproducedaninefficientsubsur- facedrainagenetwork, whichjoineda moreefficient,channelised networkdrainingwaterfromthecleanericefartherupglacier.

Consequently, while the intricacies of the englacial and sub- glacial drainage systems of Alpine and Arctic clean-ice glaciers havebeenwell-studiedanddocumented(e.g.reviewsbyHubbard andNienow,1997; FountainandWalder,1998; Irvine-Fynnetal., 2011),thereislimitedknowledgeofhowwaterflowsthroughand beneathdebris-covered glaciers.Ouraimhereinisto usefluores- centdyetracingexperimentstoinvestigatethehydrologicalsystem ofthedebris-coveredKhumbuGlacier,Nepal.Ourspecificresearch objectivesare to: i)determine whethera subsurfacehydrological systemexists within and/or beneath Khumbu Glacier;ii) if such a systemexists, determine its hydraulic characteristics, including likelyflowpaths;andiii)furtherelucidatethenatureofmeltwater transportthroughtheglacier’slinkedchainofsupraglacialponds.

2. Methods

Dyetracing experiments were carried out onKhumbu Glacier (Fig. 1) during the 2018 pre-monsoon season. The glacier has a well-documented, expanding area of supraglacial ponds that in 2015 covered 3.2% of the glacier’s 7.1 km² debris-covered area (Watson et al., 2016). Ponds are particularly prevalent along the glacier’seasternmargin,wheretheyhavebeencoalescingandcon- necting hydrologically over recent years to form a linked chain (Irvine-Fynnetal., 2017; Watsonet al., 2016).A large, perennial supraglacial channel originating in the upper clean-ice region of the glacier’s ablation area (>9 km upglacier from the terminus) has been present since at least 2005 (Gulley et al., 2009). Field observationsin2017foundthatthischannelprogressivelyincises downglacier, becoming englacial just above the confluence with Changri SharGlacier (Fig. 1). Khumbu Glacier’s hydrological sys- temreceives an additionalinput ofmeltwater fromChangriShar andChangri NupGlaciers (hereafter the Changri catchment), the proglacial streams of which coalesce and cascade down a steep gorgeinto,andlikelybeneath,KhumbuGlacier(Benn,pers.comm., 2018).AnoutburstfloodfromChangriSharGlacierin2017entered KhumbuGlacieratthispointandappearedtotransitmuchofthe glacier’s length below the surface (E.S. Miles et al., 2018). Little elseisknownabouthowwaterflowsthroughKhumbuGlacier,but currently,onlyonedominantactivesupraglacialchanneldrainsthe

glacier,formingtheproglacialstreamontheeasternmarginofthe terminus (Fig. 1B). Ananalysis ofavailable satellite imagery con- firmsthatthisconfigurationhasprevailedsince atleasttheearly 1990s.

Sixteendyetracingexperimentswereundertakenbetween27th April and 14th May 2018 across the lower ∼^{7 km of Khum-} bu’s9 km-longablationarea.Fluoresceindyewas usedduetoits photo-degradation,resultinginminimaldownstreamimpacts.Dye wasinjectedintosupraglacialstreamsorchannelsatselectedloca- tions(dye injectionpoints (DIP^∗),Fig.1),involumesof1–150 ml accordingtothe tracingdistance(Table 1). Moredyewouldhave allowed clearer breakthrough curves, however injection volumes wererestrictedtominimisedyevisibilitybeyondKhumbuGlacier duetotheglacieranditsproglacialstreambeinglocatedinaNa- tionalParknearpopulartrekkingroutes. Inall cases,dyewasde- tectedusingTurnerDesignsCyclops-7fluorometers(F^∗,Fig.1)log- gingatone-minute intervals, locatedatstrategic junctionsin the supraglacial hydrologicalnetwork. Allfluorometers were shielded fromdirectsunlightusinganinbuiltshadecap,andwerefullysub- mergedbeneaththewater.Imagesandavideo ofadyeinjection arepresentedintheSupplementaryMaterial(Fig.S1andVideoS1, respectively).

Dyetracingexperimentswerecarriedoutatfourdifferentspa- tialscalesovertheglacier:i)terminus(300 mfromtheterminus), ii) short-range (within ∼^{600 m of the terminus), iii) long-range} (∼^{7 kmfromtheterminus),andiv)pond-based}experiments(con- ductedinvarious supraglacialpondsnearthe terminus;Table1).

One terminus trace was carried out (dye injected at DIP1), with dyedetectedby asinglefluorometerontheterminusstream(F0;

Fig.1C).Threepond-basedtraceswerecarriedoutjustupglacierof thisstream(DIP2,DIP3andDIP4),alsodetectedby thesameflu- orometerwhichwassubsequentlyremoved.Sixshort-rangetraces were conducteda shortdistance upglacieragain (DIP3,DIP5 and DIP6),detectedbytwofluorometerslocatedontwodistinctinlets tothefinalsupraglacialpond(F1andF2).Thesefluorometersalso detected dyereturns fromfivelong-rangetraces(DIP7 andDIP8).

Anadditionalfluorometer(F3)was positionednearthemiddleof the linked supraglacial pond chain forthe later three long-range traces;onefinalpond-basedtracewascarriedoutabovethisfluo- rometerwithinthelinkedpondchain(DIP9).

Thefluorometerusedforthefirstsixtraces(oneterminustrace andthreepond-basedtracesatlocationF0,andthefirsttwoshort- rangetracesatlocation F1) could notbe calibratedto recordab- solute concentration.Thisfluorometerwas thereforereplacedand notusedagain;thesefirsttracesarehereafterreferredtoas‘pre- tests’ (PT^∗) due to the extra data correction that was required to set the background to ∼^{0 parts per billion (ppb; an offset of} +^{1,250 ppbwas applied). This offset did not influence any sub-} sequentanalysis. Alltestswiththeremaining three fluorometers, whichallfunctionedandwerecalibratedcorrectly,arereferredto as‘dyetraces’(DT^∗).

Measured fluorescein concentrationswere correctedfor water temperature,whichwas alsorecordedby thefluorometerloggers, asfollows:

Fr

=

.

where Fr is the calculated fluorescence atthe referencetemper- ature, Tr; Fs isthe observed fluorescenceat thetime ofreading the sample temperature, Ts;andn is thetemperaturecoefficient forfluorescein(0.0036^◦C⁻¹)(TurnerDesigns,2018).Asmallnum- ber ofextremedatapoints, includingnegative values orreadings that exceeded the maximum detection of the fluorometer, were removed andtheremaining datawereinterpolatedlinearly.Inall cases,suchoutlierscomprisedindividual,isolatedpoints,whichwe assumeresultedfromrareelectronicdisturbance.

Fig. 1.A)LocationofKhumbuGlacier,Nepal;B)thelower∼8kmofKhumbuGlacier,showingdyeinjectionpoints(DIP^∗)andfluorometerlocations(F^∗);C)thelower

∼^{500 mofKhumbuGlacier.ContoursaremarkedinB)every100 m(from4,900 ma.s.l.attheterminusto5,200ma.s.l.byDIP8).Field}observationsaremarked,including waterseepagefrom/intosupraglacialpondsbeneaththedebrislayer,a shallowmoulinandenglacialreservoirs.ThepondbelowDIP6onlyappearedfollowingfieldwork,and isinferredtohavebeenlocatedenglaciallyduringthefieldseason(seeSection4.2).SupraglacialpondsandstreamsweremappedmanuallyfromthebackgroundPlanetScope OrthoScene,capturedduringthefieldseasonon24.04.2018(PlanetTeam,2017).

Table 1

Keydatafromthebreakthroughcurvesandsubsequentanalysisofthesuccessfulpre- anddyetraces(PT^∗andDT^∗,respectively).Thelong-rangetestsreturnedtwodistinct dyebreakthroughcurves(BTC1and2).

Dyetrace type

Test

#

Dyein- jection point (DIP)

Dateand time

Dye volume^**

(ml)

Fluorometer 1^* Fluorometer 2

Distance (m)

Transit time (min)

Velocity (m s⁻¹)

Disper- sivity (m)

Distance (m)

Transit time (min)

Velocity (m s⁻¹)

Disper- sivity (m)

TERMINUS PT1 DIP1 27/04/2018 11:06 1 156 11 0.24 3.0 – – – –

POND TESTS PT2 DIP2 27/04/2018 11:37 1 234 157 0.02 – – – – –

DT7 DIP9 08/05/2018 13:48 60 450 2,510 0.003 – – – – –

SHORT-RANGE TESTS

PT5 DIP3 28/04/2018 11:05 10 82 24 0.057 3.8 – – – –

PT6 DIP5 28/04/2018 12:23 15 105 29 0.060 6.0 – – – –

DT1 DIP3 30/04/2018 10:00 10 82 19 0.072 1.3 105 31 0.056 1.2

DT2 DIP6 30/04/2018 10:43 20 193 191 0.017 11.3 212 210 0.017 14.1

DT5 DIP3 08/05/2018 11:24 20 82 38 0.036 9.1 105 49 0.036 7.8

DT6 DIP6 08/05/2018 12:06 30 193 190 0.017 15.8 212 225 0.016 20.4

LONG-RANGE TESTS(BTC1)

DT3 DIP7 30/04/2018 15:34 100 5,067 5,089 0.017 0.00004 5,077 5,097 0.017 0.0011

DT4 DIP8 06/05/2018 12:10 100 7,084 2,302 0.051 0.00037 7,094 2,310 0.051 0.0059

DT8 DIP7 09/05/2018 10:14 100 5,067 4,761 0.018 0.00005 5,077 4,768 0.018 0.0005

DT9 DIP7 14/05/2018 15:50 150 5,067 1,638 0.052 0.00038 5,077 1,646 0.051 0.0083

LONG-RANGE TESTS(BTC2)

DT3 DIP7 30/04/2018 15:34 100 5,067 5,773 0.015 0.05291 5,077 5,773 0.015 0.0645

DT4 DIP8 06/05/2018 12:10 100 7,084 2,361 0.050 1.53 7,094 2,369 0.050 0.76

DT8 DIP7 09/05/2018 10:14 100 5,067 5,501 0.015 1.25 5,077 5,500 0.015 1.24

DT9 DIP7 14/05/2018 15:50 150 5,067 2,956 0.029 0.43 5,077 2,962 0.029 0.61

NOBREAK- THROUGH

PT3 DIP3 27/04/2018 12:06 2 – – – – – – – –

PT4 DIP4 27/04/2018 12:47 5 – – – – – – – –

DT10 DIP8 15/05/2018 10:12 150 – – – – – – – –

* Fluorometer0forallpre-traces(PT^∗)andFluorometer3forDT7only.

** Dyewasa41%solution.

Fig. 2.Fullfluorometer(F^∗)timeseries(sampledeveryminute)duringthe2018pre-monsoonfieldseasonfor:A)F1;andB)F0andF2.Dyeinjectionsareshownwith orangemarkers.Darkgreybarsindicatethetimeperiodsofeachindividual,successfuldyetraceandpre-trace(PT)fromtimeofinjectiontotheendofthebreakthrough curve(BTC).Lightgreybarsshowuncertainbreakthroughs(orbreakthroughendpoints).Datawerefilteredtoremovebackgroundnoise(seeMethods),whichreducedall peakmaxima.ThebackgroundofF0andF1(PTs)hasbeenrealignedto∼^{0ppb(seeMethods)andthetimeseriesscaledforcomparisonwiththeotherfluorometerseries.}

Each dye breakthrough curve was identified and dye transit time (tm, s) calculated from the time of injection to the maxi- mum peak.Minimum transitvelocity estimates(u_m, m s⁻¹) were calculatedusingthetransitdistance(x,m):fortheterminus,pond- based,andshort-rangetraces, thestraight-linedistancewas used (Hubbard and Glasser, 2005; Seaberg etal., 1988). Forthe long- rangetraces,thestraight-linedistancewasadaptedslightlytofol- lowtheglacier’scentreline.The dispersioncoefficient(D,m²s⁻¹) wascalculated:

D

=

(

−

)

(

)

[

(

j

)

]

wheretj (s) represents t1 andt2,the time of half the peak dye concentration ontherising andfalling limbs ofthebreakthrough curve, respectively; tm (s) isthe time to peak concentration, ob- tainedby definingtheaboveequationfor j=^{1 and2andsolving} iterativelyfort_m untila commonvalue ofD isobtained(Seaberg etal.,1988; Willisetal.,1990).Dispersivity(d,m)wasthencalcu- lated:

d

=

(3) Thefluorometerdatawerefilteredusinga15 min windowfast Fourier transform to remove daytime noise, which was likely a consequence of an increased and more variable suspended sedi- mentconcentrationwithasimilarfluorescencewavelengthtoflu- orescein during hours of greater discharge (Smart and Laidlaw, 1977). No concurrent discharge data were collected to corrobo- rate this, which also precludes the calculation of dye recovery.

The rawfluorometertime series are providedin theSupplemen- tary Material (Fig. S2). Due to the limited duration of some dye breakthroughcurves,thefiltereddatawerenotsuitablefortransit velocity oranysubsequent calculations. Thefiltered curveswere, however,used toverify thevalues oftj selectedto calculatethe dispersioncoefficient,whichwereoftendifficulttodiscerndueto the backgroundnoise. The difference indispersivity betweenthe filteredandunfiltereddatawasnegligible(<1 m).

Dischargedataareavailablefromsaltdilutiontestsduringthe 2017 pre- and post-monsoon seasons (Table S1), both from the upglacier supraglacial stream(DIP8) andat the terminus(by F0;

Fig.1). Whilst acknowledging that thisdataset was collected the year before this study, exploration of the glacier’s supraglacial drainagesysteminbothyearsdidnotshowanysignificantrecon- figurationbetween2017–2018.Thesedischargemeasurementsare thereforepresentedasanapproximationoftherelativedischarges of the glacier’s supraglacial stream andoutflow atthe terminus, bothbeforeandafterthemonsoon.

3. Results

Thefull time seriesofthefluoresceinconcentrationsrecorded bythreefluorometers,F0,F1andF2,areshowninFig.2.Thirteen ofthesixteendyeexperimentsgavesuccessfulreturns,whilethe three that did not were not unexpected (Table 1): PT3 and PT4 likely merged with PT2, and the missing breakthrough of DT10 was probably caused by the fluorometers being removed before the dyeemerged. The peak maxima ofthe long-rangetraces are substantiallysmaller(>6×^{smaller)thanthoseofthe}short-range andpond-basedtraces.Thisisunsurprisinggiventheconsiderable distances travelled forthe long-range traces (∼⁷⁰× ^{greater than} theshort-range traces) whiledyevolumeswere only scaledby a factorof∼^10.

A single terminus trace was conducted to obtain a measure- ment of the outlet stream velocity and dispersivity (Fig. 3). The

Fig. 3.Fluorometer(F0)timeseriesfromtheterminuspre-trace(PT1;dyeinjection point(DIP1)).ThebackgroundofF0hasbeenrealignedto∼0ppb(seeMethods).

breakthrough curve shows a short transit time with a narrow breakthrough curve width,implying low dispersivity. Thisissup- portedbythetestdata(Table1),indicatingavelocityof0.24m s⁻¹ andadispersivityof3.0m.Thevelocitiesoftheshort-rangeexper- iments are an orderofmagnitude lower thanthe terminustrace (Table 1; Fig. 4), but with similarly low dispersivity values. The breakthroughs recordedatbothfluorometersfromDIP3andDIP5 tend tobe fasterandlessdispersed thanthose fromDIP6.These traces alsoshow adifference in peaktiming betweenF1 andF2, with breakthrough maxima at F2 occurring 11–35 min after F1.

The final short-range traces (Fig. 4C) show alower velocity from DIP3,anda muchhigherdispersivityforboth injections(∼^{8.0 m} increaseforDIP3,and∼^{6.0 mincreaseforDIP6relativetothepre-} vioustrace;Table1).

Thefluorometerreturnsfromthepond-basedtracesareshown in Fig. 5. All four breakthroughs are uncertain andPT3 and PT4 areindeterminablehavingmergedwithPT2(Fig.5A).Wherepeaks seemapparent (particularlyforDT7;Fig. 5B), theverysmallcon- centration rangeshould be considered alongside theshort trans- portdistancesforthesetraces(andthebackgroundstreamfluores- cence,alsoshowninFig.S2).Velocitieswerecalculatedtothefirst major peak ofPT2 (0.02 m s⁻¹) and DT7 (0.003 m s⁻¹; Table 1), butdispersivitiescouldnotbecalculatedduetotheindistinctend points ofallbreakthroughcurves.Fluorometer3wasleft running until the endofthe measurement periodbutshowed onlyback- groundfluorescencevalues(seeFig.S2).

The breakthrough curves for the long-range experiments are presented in Fig. 6. All traces from DIP7 returned two dis- tinct breakthroughs recorded at both fluorometers, separated by 10–20 h (Table 1; Fig. 6A, C, D). Two breakthroughs may also have been returned from the DIP8 trace (Fig. 6B; clearer from the raw data) but the period between the breakthroughs was much shorter (<1 h). The first breakthrough consistently pro- duceda peakofshorterduration, whilethesecond breakthrough was slower (by 0.001–0.003 m s⁻¹ over 5–7 km) and more dis- persed.Averagetracevelocities(ofbothbreakthroughs)aresimilar to the short-range traces(∼^{0.03 m s−1). Thedifference in break-} throughs betweenthe fluorometers was present butlessdistinct compared to the short-range experiments: the F2 breakthrough maximawere,onaverage,5.5 min laterthanforF1.ThefinalDIP7 trace(Fig.6D)occurredinapproximatelyhalfthetimeofthefirst twotraces.

Fig. 4.Fluorometer(F^∗)timeseriesfor:A)theshort-rangepre-traces(PT^∗);and(B, C)theshort-rangedyetraces(DT^∗).Thedyeinjectionpoint(DIP^∗)andtimeisnoted foreachexperiment.ThebackgroundofF1(PTs)hasbeenrealignedto∼0 ppb(see Methods).

SummaryvelocityanddispersivitydataarepresentedinFig.7.

Thedifferencesbetweentheextremesoftransitdistancearehigh- lighted at either end of the plot: the shortest distance trace (DIP1) has the fastest velocity and a low dispersivity (156 m, 0.24 m s⁻¹ and3.0 m, respectively; Table 1); the longest transit distance(DIP8) showslow velocities anddispersivities (7,084 m;

<0.018 m s⁻¹ and<1.25m,respectively). Therepeat short-range experimentsshow amoderatedecreaseinvelocityandsignificant increaseindispersivity overtime (e.g.from0.072 to0.036 m s⁻¹ andfrom1.3to9.1 mforDIP3,respectively).Therepeatlong-range tracesfromDIP7show athree-fold increaseinvelocity (e.g.from 0.017 to 0.052 m s⁻¹ forthe first breakthroughof DT3 andDT9, respectively).

4. Discussion 4.1.Terminustrace

Thesingle terminusexperiment was conductedto aid thein- terpretation of all subsequent dye traces. Velocities ≥^{0.2 m s−1} and dispersivities ≤^{10 m have} traditionally been interpreted to show drainage through an efficient, fast-flowing, integrated and channelised drainage system (Burkimsher, 1983; Seaberg et al.,

Fig. 5.Fluorometer(F^∗)timeseriesfromthepond-basedexperiments:A)pre-traces (PT^∗)towards the terminus;and B) dyetrace(DT7)withinthe linked chainof supraglacialponds.Thedyeinjectionpoint(DIP^∗)andtimeisnotedforeachex- periment.Datawerefiltered(boldlines;seeMethods)toremovebackgroundnoise.

Raw,unfiltereddataareshownwith athinner,moretransparentline.Theback- groundofF0hasbeenrealignedto∼^{0 ppb(seeMethods).Thefallinglimbatthe} startoftheF0timeseriesinA)isthatofPT1(cf.Fig.3).

1988; Williset al., 2012, 1990).Both themoderately highveloc- ity (0.24 m s⁻¹) and low dispersivity (3.0 m) of this test con- firm our observations that this large, single supraglacial channel evacuates meltwater rapidly froma smallsupraglacial pond into the proglacial stream (approximate discharge 1–2 m³s⁻¹, Table S1; image of injectionandchannel shownin Fig. S1A). This dye trace therefore provides a reference for the maximum efficiency wemightexpectwithinthesystem(Fig.7).

4.2. Short-rangetraces

As the only visible exit for meltwater on Khumbu Glacier is fromthesupraglacialsystematthe terminus,theshort-range ex- periments were carried out to characterise this sector, which is alsothe endof thelong-rangedrainage system. Thevelocities of 0.02–0.07 m s⁻¹ indicateslowtransit,withconsistentlyfasterand lessdispersedtransportfromDIP3thanfromDIP6.Directfieldob- servations indicatedthat DIP3was locatedattheheadofa short lengthofsupraglacialchanneldrainingasupraglacialpond,which flowed into a shallow moulin feeding a larger, shallowenglacial reservoir(Fig.1C;Fig.S1C–D).Thisreservoirwasvisiblethrougha surfacefracture,flowingintothesmallsupraglacialpondaboveF1.

Afterinjection,dyecouldbefolloweddownthestreamandviewed intheenglacialreservoir throughasurface fracture,only slowing in the supraglacial pond immediatelyabove F1.To produce slow velocitiesbutreasonablylowdispersivities,weinferthatthereser- voirandsupraglacialpondprovidedtemporarystoragebuttransit was relatively undisturbed by flow complexities, such as turbu- lence and/or eddies,limiting dispersion. In addition to the main outflowpast F1,dyewas observedtoleave thispond by seepage underthedebrisalongtheeasternpondmargintowardsthenext downstreampond(Fig.1C),aprocessthathaspreviouslybeensug- gested but not observed (Irvine-Fynn et al., 2017). The dye was

Fig. 6.Fluorometer(F^∗)timeseriesfromthelong-rangedyetraces(DT^∗)fromDIP7(A,CandD)andfromDIP8(B).Thedyeinjectionpoint(DIP^∗)andtimeisnotedforeach experiment.Datawerefiltered(boldlines;seeMethods)toremovebackgroundnoise.Raw,unfiltereddataareshownwithathinner,moretransparentline.Therawdata peakinpanelDhasbeencut-offtoharmonisethey-axisscales.

Fig. 7.Summaryvelocity(bluesymbols)anddispersivity(redsymbols)datafrom eachtesttype(cf.Table1).Eachdyeinjectionpoint (DIP^∗)hasaseparatesym- bol:forthelong-rangetests,thelargersymbolindicatesthefirstdyebreakthrough, andthesmallersymboltheseconddyebreakthrough.Wherethedyebreakthrough wasrecordedbytwofluorometers,valueswereaveragedtoshowthe difference betweenthedyebreakthroughsmoreclearly.Thedottedlinesandshadedareasin- dicatevaluesthatdeterminewhenadrainagesystemistraditionallyconsideredto bechannelisedratherthandistributed(velocity≥^{0.2 m s−1and}dispersivity≤^{10 m,} Burkimsher,1983; Seabergetal.,1988; Willisetal.,1990).

well diffusedby theslow, constant flowthrough thesupraglacial pond; it isassumedthat our breakthroughs fromthemain pond outlet, on which F1 was located, are representative of all pond outflows. The multiple visible and concealed outflows from this one small pond highlight the distinctive near-surface complexity ofdebris-coveredglacierhydrology.

TheslowervelocitiesandgreaterdispersivitiesfromDIP6com- pared to DIP3 imply greater storage in the 100 m upstream of

DIP3(transittimeofhoursratherthanminutes).DIP6waslocated at the downglacier end of a large supraglacial pond, at the end of the linked pond chain. After injection, dye moved slowly to- wards the down-flow tip of the pond, exiting the pond beneath the debris (image in Fig. S1B). How water was transported be- tweenDIP6andDIP5isnotknown,butweinferthatthisflowpath included an additional, larger, englacial reservoir(s) dueto: i)an absence ofsurface waterin thevicinity; ii) thefar slowerspeed andgreater dispersionoftheDIP6 breakthroughcurverelative to DIP5andDIP3;andiii)theonlyvisibleinputtothepondatDIP5 beingasmallinlet,alsobeneaththedebrislayer.Satelliteimagery showstheformationofasupraglacialpondbetweenDIP6andDIP5 shortlyfollowingthe fieldexperiments(Fig. 1C),whichmayhave resulted fromtheflooding and/orsurface collapseofan englacial reservoir.

Thedifferenceinthetimeofbreakthroughatthetwofluorom- eterswasconsistentforalltheshort-rangeexperiments,revealing a slight variation indrainage betweenF1 andF2 downstream of DIP3. The drainage network leading to F1 was noted above, but the systemfeeding F2 is unknown,despite exploration. We infer thataportionoftheflowisdivertedfromtheF1pathwithinthe englacialreservoir,andfollowsalessdirectroutebeforeemerging as theF2 stream. Due to thedispersivity variations between the short-range tests, weare unable toconfirmwhetherthe network comprisesanadditionalenglacialreservoir,asinuousconduitora smalldistributedchannelnetwork.

The repeatshort-range tracesrevealadecrease insystemeffi- ciency over thefield season, withthe final experimentrecording 50% of the velocity measured in previous tests, and dispersivity values 140–700% greater (Table 1; Fig. 7). Although the limited number oftracesshould beacknowledged, theseresultsdo align withchangesobservedinthefield.Beforethefinalrepeattraceon 8thMay,wefoundthatthemoulinbelowDIP3hadcollapsedand re-routed flow intoa shallowergradientmoulin ∼^{5 m upstream.}

Inthissection,thestreamhadagreaterdischarge andfasterflow than had been previously observed. As the dye still reached F2, the system changes didnot cut thisdrainage route off, confirm-

ing that the flow divergence was downstream of thesemoulins.

Thebreakthroughcurves,particularlyfromDIP3,display multiple peaks(Fig.4),suggestingthatthisnewpathwayroutedwaterless efficiently through multiple conduits into the englacial reservoir, despitethe greater discharge intothe moulin.Debris inthecon- duit may have contributed to the higher dispersivity (Gulley et al., 2014), butthe newly-formed pathwaymaysimply havebeen moreconvoluted,possiblyduetoitexploitingvoidsorpre-existing weaknessesinthe ice. Dye mayalsohavebeen delayedanddis- persed in the englacial reservoir due to its higher entry point.

Ourobservationsofthesepathwaysandtheirrapidchangesinhy- draulictransfer displaysa further hydrologicfeature that maybe distinctiveto debris-covered glaciers: re-routing occurringdueto surface collapse, which is extremely prevalent on debris-covered glaciersduetotheirhighly spatiallyvariableratesofsurface low- ering(Bennetal.,2017;E.S. Milesetal.,2017a,2018).

4.3.Pond-basedtraces

Thedyeforthepond-based experimentsattheterminus(PT2 toPT4;Fig.5A)wasinjectedover70 min betweenthethreeDIPs, but there is only one clear breakthrough maximum (from PT2), yielding a low transit velocity (0.02 m s⁻¹). The half-life of fluo- resceindepends ona large numberoffactors, butis shortwhen exposed directly to bright sunlight (Smart and Laidlaw, 1977):

smallsample tracessuggest completephoto-degradationwithin a few hours (Turner Designs, 2018). When dispersed in more tur- bidwater,through whichlight isattenuated, the half-lifewill be correspondingly longer. We suggest that dyetransit was so slow through the ponds that the separate dye injections merged and muchofthedyephoto-degradedinthelargeice-freepondimme- diatelybeyondF1/F2.Thesameconclusioncanbedrawnfromthe otherpond-basedtest(DT7fromDIP9;Fig.5B),theindeterminate breakthroughofwhichproduced avelocity of0.003 m s⁻¹ tothe largestfirstpeak, ∼^{40 h afterinjection(cf.Fig.S2).Thisisanor-} derofmagnitudeslowerthanterminuspondtracesandtwoorders ofmagnitudeslowerthantheterminustrace(Fig.3).The ratioof throughflowingdischargeto pond volumeistherefore inferredto beexceptionallyslowpriortothemonsoonmeltseason,support- ing previous findings (Irvine-Fynn et al., 2017). No dye was de- tectedatF1orF2;itmayhavebecomeverywelldispersedwithin thepondand/or beenlargelydestroyedby photo-degradationbe- tweenF3andDIP4.

4.4.Long-rangetraces

Although of low concentration, the long-range tests yielded breakthroughcurves (Fig. 6), confirming that meltwateris trans- portedfromthe upperablation area to the terminusof Khumbu Glacier,exitingtheglacieratthesurfaceratherthan(entirely)be- inglosttogroundwater.Theslowvelocitiesforthelong-rangetests (mean of ∼^{0.03 m s−1) and} unexpectedly low associated disper- sivities (all <1.5 m; Table 1, Fig. 7) indicate minimal long-term storageand/oreddyingwithinthetransportpath,suchasmightbe foundinsupraglacialorenglacialponds.Theabsenceoflong-range trace breakthroughs at F3 further indicates that this subsurface drainagesystemdoesnotlinkintothesupraglacialpondchain.

Compared to previous dye tracing studies on other debris- covered and clean-ice glaciers that inferred flow through an ef- ficient, channelised drainage network partly from fast through- flowvelocities (Burkimsher,1983; Fountain,1993; Hasnain etal., 2001; Nienow et al., 1998; Pottakkal et al., 2014; Schuler et al., 2004; Seaberg et al., 1988), our trace velocities are one, and in somecasestwo,orders ofmagnitudelower.Ourresultsaremore akin, for example, to those from the similar-sized but debris- free Midtdalsbreen, Norway, which were interpreted in terms of

drainagethroughalinked-cavitysystem(Willisetal.,1990).How- ever,thesenetworksproducedmarkedlyhigherdispersivityvalues (10.0–71.1 m)than werecord atKhumbuGlacier.Itthereforeap- pears that at least part of the subsurface drainage network at KhumbuGlacierdiffersfromeachofthetwoformspreviously re- portedfromclean-iceglaciers:inshowingarelativelylowvelocity andlow dispersivity itis neithercompletely channelised/efficient (highvelocityandlowdispersivity)norcompletelydistributed/in- efficient(lowvelocityandhighdispersivity).

Wepropose thatKhumbu Glacier’sdrainagesystemacross the ablation area initiates supraglacially in the clean-ice region be- neath the ice fall above DIP8. Gulley et al. (2009) mapped sev- eral supraglacialchannels incisinginto ‘cut-and-closure’ englacial conduits in thisregion in2005, thepathways of whichare very similar to parts of the sinuous supraglacial channel network we observedin2017–18(Fig.1B; imageinFig.S1F). Wesuggestthat rapidsurface melthascausedtheseenglacialconduitstobecome re-exposed atthesurface, ashasbeen observedon other debris- covered glaciers (Miles et al., 2017a). Indeed,this allowed usto hike along the progressive incision ofthe supraglacial stream in 2017 to where the stream now disappears to become englacial,

∼^{300 mupglacier ofthelargemeltwater inputfromthe Changri} catchment. At this point, multiple relict conduits were visible above the active channel, supporting the interpretation of an in- cisingcut-and-closureconduit.Thelocationofthestreamsubmer- gencetobecomeenglacial,andthesimilaritiesofthetracebreak- throughs fromDIP7 and DIP8,lead usto infer that this channel likelycontinuestoincisedownwardstojointhestreaminputfrom the Changricatchment. TheChangri inputwas observed toreach andfollowthebedofKhumbuGlacierin2006(Benn,pers.comm., 2018), and we expect that this is still the case due to the sub- stantialdischargeofthistorrentduringourfieldobservationsboth beforeandafterthemonsoonseason(Fig.S1Eshowsanimage of dye injectionatDIP7). This hasallowed the subglacial streamto adopt a stable position at the bed in this section of the glacier (Benn,pers.comm.,2018).

KhumbuGlacier,likemostglaciersinthisregion,hasalargely impermeable terminalmorainethat haspreviouslybeen notedto providea highlocalbase-levelforthe glacier’senglacialdrainage (Gulley et al., 2009). On the basis of the dye breakthroughs at the surface near theterminus and onlyone supraglacialchannel draining the glacier, we suggest that the highhydrological base- level for the lower glacier has prevented the subglacial channel continuingalong thebedfurtherdownglacier,resultinginan‘up- routing’ of the system back to the surface (of no more than a few tens of metres, Gades et al., 2000). This would produce an englacial depth limit to the cut-and-closure mechanism in the lower part of the ablation area – possibly following the cold- temperate ice boundary inferred by K.E. Miles et al. (2018). Be- tween thesubglacial–englacial transitionandthereappearance at the surface, the drainage continues beneath the surface, bypass- ing thesupraglacial andnear-surfacehydrological networkinthe lowerablationarea.

The subsurface drainage system downglacier of the Changri Shar confluence may therefore be a continuation of Khumbu’s sinuous,upglaciercut-and-closurechannel.Speleologicalinvestiga- tionsontheneighbouringNgozumpaGlacierhaveshownthatcut- and-closure englacial conduits can produce highly sinuouschan- nels (Bennetal.,2017; Gulleyetal., 2009). Indeed,thismayalso be prevalent within Khumbu Glacier,due to the glacier’sdebris- coveredtonguehavingaverylowsurfacegradient(Fig.1B).These may combine to produce a low hydraulic gradient and thus en- courage meandering.Asinuoussystemcould potentially resultin the actual transit velocities being up to an order of magnitude higher (using theactual, rather thanstraight-line distance): such velocitiesmayevenbesimilarto thoseusedtoinferchannelised,

ratherthandistributed,drainageonclean-iceglaciers(e.g.Nienow etal.,1998).Theywouldalsocorrespondbettertotheassociated tracedispersivitiesindicatingundisperseddyetransitthroughthe drainage system (Fig. 7). However, such velocities would still be lower than those observed on the debris-covered Miage Glacier, which hasa much steeper ablation area, potentially encouraging faster flow (Fyffe et al., 2019), and the overall transit through Khumbu Glacier is still slowand inefficient. We suggest there is a short section beneath the Changri Shar confluence that is not sinuous(Benn,pers. comm.,2018),havingadapted totheChangri catchmentinput.Anotable,straight,surfacedepressionfollowsthe western glacier margin for a kilometre or so, beginning imme- diately belowthe confluence,which maybe similar to the sub- marginal channels observed on Ngozumpa Glacier (Benn et al., 2017; Thompsonetal.,2016).

The dual breakthroughs for our long-range traces are most likely caused by a division of the drainage system (Burkimsher, 1983; Willisetal.,1990).Giventhatweinterpretthisflowaspre- dominantlyenglacial,thesedualpathwayscouldbesimilartothe englacialsidepassagewaysreportedonNgozumpaGlacierbyBenn etal. (2017).

Further complexity in Khumbu’s subsurface drainage system is suggested by the 2017 pre-monsoon salt dilution discharge measurements (Methods; Table S1). At 13:00on 24thMay 2017, the discharge ofthe supraglacial stream nearDIP8 was ∼^1.5m3 s⁻¹; twodayslater, alsoat13:00,the discharge oftheproglacial streamnearF0was ∼^{1.3 m3s−1.Comparable}measurementswere found in the post-monsoon season: on 22^nd October 2017, the supraglacial stream discharge was ∼^{2.4 m3s−1; two days later} at 13:00, the proglacial stream discharge was ∼^{2.6 m3s−1. The} similarityofthesevalues indicates thatthe principalsupraglacial streamalonecouldcontribute almostall oftheglacier’sproglacial discharge,both pre- andpost-monsoon.Khumbu’s fullsubsurface meltwaterdischarge,includinginputsfromtheChangricatchment, maynot thereforeemerge fromtheglacier’sterminal portal; yet, nonotableupwellingemergesfromtheglacier’sterminalmoraine.

Itispossiblethatthe2017dischargemeasurementsdonotreflect longer-term conditions.Alternatively, ifthese discharge measure- ments are representative more generally, then a sizeable propor- tion of Khumbu’s subsurface meltwater must emerge elsewhere fromthesystem.Wenotethatthiscannotholdforalltheglacier’s subsurfaceflowbecausethelong-rangedyetracesgavesuccessful breakthroughs atthesurface neartheterminus. Thelost compo- nentcould bestoredwithin theglacier,flowasgroundwaterthat emergesfartherdown-valleythantheterminalmoraine,oremerge diffusely acrosstheouter slopeofthe terminalmoraineata rate thatisinsufficienttoovercomelocalevaporation.Allofthesepro- cesses deservefurtherinvestigation.The discharge measurements indicatefurtherthatthelinkedsupraglacialpondchainprovidesa relativelysmallproportionoftheglacier’sdischarge.

The bypassing of the supraglacial hydrological system by the subsurface system for much of the lower ablation area is simi- larto the pathwayobserved for perennial sub-marginal conduits on Ngozumpa Glacier, which alsoroute back to the surface very close to the terminus into the proglacial Spillway Lake (Benn et al., 2017). Khumbu Glacier’s subsurfacedrainage isroutedto the surface betweentheendof thelinked pond chain(F3)andF1/2.

Thetimedifferencebetweenthefluorometerbreakthroughsissim- ilarto thatoftheshort-range traces(F2 beingconsistently later), confirming that the flow into F1/2 only diverges very close to the fluorometers and that the subsurface drainage system joins theshort-rangedrainagesystemweobserve,perhapsnear/intothe englacialreservoirvisiblefromthesurface.However,theupglacier subsurfacedrainagesystemmaynotbeperpetuallyseparate from the supraglacial system: evidence for the periodic drainage of perchedponds into the englacial network hasbeen observed on

both KhumbuandNgozumpaGlaciers,aswell asatother debris- covered glaciers (Benn et al., 2017, 2009; Miles et al., 2017a).

Further, E.S. Miles et al. (2018) observed a lake outburst event that propagated through Khumbu Glacier in 2017, suggestingan overflow pathwaycreated by theflood waters emerged fromthe subsurfacedrainagesystemintothelowersupraglacialpondchain.

The breakthroughforthefinal tracefromDIP7(DT9)occurred notably sooner thanthefirst twotraces(DT3andDT8), withthe firstdyebreakthroughshowingavelocityoverthreetimesgreater (Fig. 6D; Fig.7). On thenights of the1st and8thMay (5 and6 daysbeforeDT4 andDT9,respectively),therewere largesnowfall events.WethereforeinterpretthatthefastervelocitiesforDT4and DT9mayhavebeenashort-termsystemresponsetogreatermelt- waterinputs.Snowmeltfromthesecondeventbeganaroundmid- dayonthe9thMayontheglacieritself;thesnowcovermayhave reduced surface ablation,resulting in the slower velocity of DT8 thanDT4.Theinfluenceonthesubsurfacedrainagesystemwould havebeengreaterseveraldayslaterduetothedelayedmeltingof thesnowathigherelevations,contributingtothegreater velocity ofDT9.

5. Summaryandconclusions

To our knowledge, this study reports the first successful dye tracing experiments at a debris-covered glacier in the Nepal Hi- malaya andthefirstdyetracing-basedinvestigationexploring the intricacies of debris-covered glacier hydrology, both at and be- neath thesurface. Weconductedsixteendyetracing experiments on KhumbuGlacier,Nepal,whichreveal previouslyunknown fea- turesandcomplexitiesinthesurfaceandsubsurfacehydrologyofa debris-covered glacier,manyofwhichdifferfromcurrentnotions of clean-ice glacier drainage. We highlight the following conclu- sions:

• ^{A likely highly sinuous and} channelised subsurface drainage system exists at Khumbu Glacier, flowing for some distance alongtheglacier’sbedbelowtheconfluencewithChangriShar Glacier.Thesystemdoesnotappeartoinvolvelong-termstor- age, and may have the potential to transport water rapidly, particularlyafterheavyprecipitationevents.

• ^{Flowthroughthelinkedchainof}supraglacialpondsalongthe easternmargin ofKhumbuGlacierisextremelyslow(velocity

∼^{0.003 m s−1)andcomprisesonlyasmallproportionoftotal} flowfromtheglacier,duringthepre-monsoonseason.

• ^{Subsurface flowis a parallelsystemthat bypasses thelinked} supraglacialpondchain.Itisultimatelyre-routedback tothe surface close tothe terminus,where itjoinsthe supraglacial system feeding the proglacial stream. A proportion of the meltwaterinputstothe subsurfaceflowdonot appeartore- emerge by this route, but the ultimate destination remains unclear.

• ^{Weobservepathwaychangesinthe}short-rangelinkedsupra- glacial andshallow englacialdrainageroutes nearthe termi- nus,triggeredbychannelcollapseinducedbycontinuedtopo- graphicevolutionofthesurface.

Whilenotingthatthenumberandlengthofthetraceswaslim- itedduetothedifficultiesofworkinginsuchanenvironmentand thecomplexity ofthesystem, wesuggestthere isgreatscopefor futureinvestigationsofboththesurface andsubsurfacehydrolog- icalsystemsofKhumbu Glacierandother debris-coveredsystems inthegreaterHimalayaandbeyond.Ourresultsareinfluencedby the timing of our experiments early in the melt season: much more could be learned by repeating these tests during or after themonsoonseason whenthesubsurfacehydrology maybecome moredevelopedaftersustainedlargeinputstothesystem.