Introduction
The simulation of jökulhlaups, also known as glacial lake outburst floods, has become a standard test case for channel-flow models (e.g. Reference NyeNye, 1976; Reference Spring and HutterSpring and Hutter, 1982; Reference ClarkeClarke, 2003). Empirical results for validating the model outputs are often sparse, consisting, for example, of the discharge data from the proglacial stream and the lake; sometimes the lake and outlet temperature are also measured. However, these measurements are not sufficient to constrain all the free parameters of jökulhlaup models. Recently, efforts have been made to collect more accurate and comprehensive data on glacial lake outburst floods: Reference AndersonAnderson and others (2003) on Kennicott Glacier, Alaska, USA, and this study on Gornergletscher, Switzerland (Reference Huss, Bauder, Werder, Funk and HockHuss and others, 2007; Reference Sugiyama, Bauder, Weiss and FunkSugiyama and others, 2007, Reference Sugiyama, Bauder, Huss, Riesen and Funk2008; Reference Walter, Deichmann and FunkWalter and others, 2008, Reference Walter, Clinton, Deichmann, Dreger, Minson and Funk2009).
The data we present were collected in 2006, when the ice marginal lake Gornersee drained via supraglacial overspill. However, the first 1.5 days of the drainage were very similar to an onset of a subglacial outburst flood: the drainage moulin was filled to its top because its capacity was not sufficient to carry all the supplied lake water. Thus pressurized-flow conditions prevailed in the entire lake outlet channel during these 1.5 days. This onset period is the focus of the present investigation, during which we measured the flow speed of the lake water by means of direct tracer injections into the lake outlet. Thus, for the first time, it is possible to compare water flow speeds simulated with a jökulhlaup model to direct measurements.
We present measurements of tracer transit speeds, sub-glacial water pressure head, lake and proglacial discharge, lake temperature and meteorological conditions. From the measured transit speeds we infer water flow speeds in the connection channel between the lake and the main subglacial drainage channel. These inferred speeds are compared to modelled speeds calculated using Reference ClarkeClarke’s (2003) jökulhlaup model which is based on the Spring–Hutter equations (Reference Spring and HutterSpring and Hutter, 1982). The model is driven by the measured water pressure at both ends of the channel connecting the lake to the main drainage channel and calibrated against lake discharge and inferred water flow speeds. To fit the model to the inferred flow speeds, simulations are performed using a sinuous channel or increased heat transfer.
In the accompanying (Part I) paper (Reference Werder, Loye and FunkWerder and others, 2009), the focus lies on the influence of the jöklulhaups in 2005 and 2007 on the glacial drainage system as a whole, also using tracer experiments as the main investigative tool. The outbursts in those two years proceeded subglacially and thus were shorter and more intense. In those years, we were not able to trace the lake water directly. Instead, moulins down-glacier of the lake were used for injections, giving us the opportunity to study the overall reaction of the glacial drainage system to a large perturbation.
Terminology
We now elucidate a few concepts and define terms to allow us to describe and discuss the tracer experiments and the accompanying model and to compare them. We assume that tracer and water travel at the same velocity and, thus, the following definitions apply to both. We include a generic variable name, if it is used later. The definitions here are a slight generalization of the ones used in the companion paper (Reference Werder, Loye and FunkWerder and others, 2009).
The flow path is the path traversed by the tracer with an associated flow path length, l. The residence time, t, between two locations on the flow path is the time interval between the passage of the maximal tracer concentration at those locations. The flow speed, v, between those two locations is v = l/t, i.e. the flow path cross-section averaged speed. The transit distance is the shortest possible horizontal travel distance of the tracer between two locations and has an associated transit path (see fig. 1 in companion paper: Reference Werder, Loye and FunkWerder and others, 2009). The ratio of the residence time and transit distance gives the corresponding transit speed, 
.
Fig. 1. Overview map of Gornergletscher. Contours of the hydraulic potential, ϕ = z b + f float ρ i /ρ w(z s − z b), are plotted, where z b and z s are the bed and surface elevation, respectively, ρ i and ρ w are the densities of ice and water and f float is the fraction of the water to ice overburden pressure, set to 0.9 (Reference Flowers and ClarkeFlowers and Clarke, 1999). z s was obtained from photogrammetry and z b from radar measurements (Reference Sugiyama, Bauder, Huss, Riesen and FunkSugiyama and others, 2008). The black arrows indicate a likely location for the main drainage channel. The grey arrow shows a possible location for the connection channel from the lake to the main channel, and the grey shaded area indicates the range of possible connection channel locations considered. Moulins (triangles), boreholes (squares), lake border (solid curve), glacier border (dashed curve), maximal lake extent (dashed curve) and active seismic shots (crosses) are also marked.
This careful distinction is necessary, as with the presented field experiments only the residence time, transit distance and speed can be determined, whereas theories of glacial hydraulics work with the flow path length and speed. Only the residence time is applicable to both, so care must be taken when comparing measurements to model results; in general, the transit distance is shorter than the flow path length, due to the geometry and sinuosity of flow path and to the vertical distance covered, meaning that the transit speed is a lower bound on the average flow speed.
Field Site and Methods
Intensive field measurements on Gornergletscher were conducted in the years 2004–08, to investigate the jökulhlaups of Gornersee. Gornergletscher is the second largest glacier in the Alps (∼60 km2, 4600–2200 m a.s.l., 14 km long). Gornersee is an ice marginal lake situated in the confluence area of the two main tributaries of Gornergletscher (Fig. 1). It has a volume of ∼4 × 106 m3 and usually drains as a subglacial outburst flood over the course of 2–7 days at the beginning of summer. The lake is located at 2530 m a.s.l. and lies 6.5 km up-glacier from the terminus. The maximal ice thickness of 450 m is found 1 km down-glacier of the lake. The field site is described in more detail by Reference Huss, Bauder, Werder, Funk and HockHuss and others (2007) and Reference Werder, Loye and FunkWerder and others (2009).
Field methods and data processing
Fluorescent dyes Uranine and Rhodamine WT were used for the injections. The detection was performed at the gauging station of Grande Dixence SA, 1.25 km below the terminus of the glacier where proglacial discharge was measured. The dye was monitored continuously with both a flow-through and a submersible fluorometer. We present a total of 11 dye injections which were conducted into the lake-drainage moulin in the morning (0830 h local time, UTC+2), early afternoon (1400 h) and evening (2100 h) between 28 June and 18 July 2006. The diurnal variability of tracer experiments conducted using the same moulin can be large (Reference Schuler, Fischer and GudmundssonSchuler and others, 2004); to minimize the influence of such effects, only experiments conducted at the same time of day were compared directly.
We present the tracer transit speed; to characterize the breakthrough curves, an advection–dispersion model with storage (ADSM) was fitted to them (Reference Toride, Leij and van GenuchtenToride and others, 1999). (See Reference Schuler, Fischer and GudmundssonSchuler and others, 2004, for a glaciological application of this model.) The ADSM takes as input the concentration time series and the transit distance. The ADSM returns estimates of the mean tracer transit speed (not presented), the dispersion, D, the fraction of mobile water, β, and the exchange coefficient between mobile and immobile water (not presented). All the breakthrough curves considered in this paper were fitted accurately by the ADSM, and are thus fully described by these parameters. Furthermore, the fraction of returned tracer mass, M, is obtained by integrating the tracer concentration multiplied by the proglacial discharge and divided by the injected tracer mass.
Here we are interested in the time the tracer spends within the glacier, the glacial residence time. Hence it is important to estimate the proglacial residence time, i.e. the time spent in the proglacial stream. Eighteen injections were conducted in the proglacial stream at different discharge levels throughout the field campaign. To the resulting proglacial transit speeds we fitted the expression:
with a least-squares regression. The parameters were estimated to be k = 0.66 m1−3e s e −1, Q 0 = 6.4 m3 s−1, e = 0.34 with a coefficient of determination, R 2, of 0.8. The zero proglacial transit speed limit as Q → 6.4 m3 s−1 is unphysical, but the proglacial discharges relevant here are all above Q = 10 m3 s−1. From the fit of Equation (1) to the data, we estimated the error in the calculated proglacial residence time to be ∼30%. The measurement error on the (total) residence time was <1%. This led to a 4% error in the calculated glacial residence time, because the proglacial residence time was around 10 times shorter than the total residence time. The ADSM parameters were not corrected for possible influences of the proglacial stream.
Subglacial water pressure data were obtained from pressure transducers (Gekon 4500) installed in boreholes drilled to the bed (labelled BH1–BH4 and BH6 in Fig. 1). The lake level was monitored with a pressure transducer (Keller DX series). Lake discharge was determined from measurements of the channel cross-section and of the water flow speed (with a current meter Schiltknecht MiniWater20) and the error was estimated to be 10%. The measured discharge agreed well with the independent estimate derived from lake level, bathymetry and water input into the lake calculated with a melt model (using the method of Reference Huss, Bauder, Werder, Funk and HockHuss and others, 2007). A thermistor was immersed into the lake water directly at M4 to measure the temperature. The air temperature and precipitation were measured at an automatic weather station located off-ice at the northern glacier margin (see fig. 1 of Reference Werder, Loye and FunkWerder and others, 2009). To determine the geometry of the moulin shaft of M4, active source seismology was conducted 2 weeks after the presented tracer experiments by lowering explosives into M4. The detonations were detected on the seismic network (Reference Walter, Deichmann and FunkWalter and others, 2008) installed in the vicinity, on the glacier (personal communication from F. Walter, 2008).
Observations
In 2006, Gornersee filled until its shore reached the moulin M4 (Fig. 1). Prior to the lake drainage, M4 had a diameter of 0.5 m and carried a maximal discharge of <0.1 m3 s−1. Once the lake shore reached M4, the lake started draining into the moulin. The water level in M4 was equal to the lake level during the first 1.5 days of the lake drainage (henceforth called the onset period), which lasted from 5 July 1400 h to 6 July 2400 h. During the onset period, lake discharge into the moulin increased from ∼0.1 to 3.5 m3 s−1 (Fig. 2e). At the end of the onset period, M4 had adjusted its capacity and its water level dropped, so the lake discharge became limited by the height of the spillway. This terminated the onset period, the time during which the lake outlet was fully pressurized. Subsequently, the lake lowered its level by incising a canyon into the ice (Reference Raymond and NolanRaymond and Nolan, 2000), it subsided slightly more than 1 m d−1 and lake discharge was in the range 2–5 m3 s−1. It took about 3 weeks to empty the lake, which initially contained ∼4 × 106 m3 of water. At the end of the lake drainage, the canyon was ∼200 m long, ∼5 m wide and up to 50 m deep; the diameter of M4 reached ∼10 m.
Fig. 2. (a) Glacial transit speed, 
; (b) ADSM parameters dispersion, D; (c) the fraction of mobile water, β; and (d) the fraction of recovered tracer mass, M, from injections done at 1400 h into M4. (e) Lake discharge into M4. The two vertical lines delimit the onset period (of the lake drainage, 5 July 1400 h to 6 July 2400 h), the 1.5 days during which pressurized channel flow conditions prevailed in the moulin.
In contrast to a subglacial lake drainage, the lake outlet (M4) was accessible and we used it for tracer injections and direct discharge measurements; the results are shown in Figure 2. The glacial transit speed was steady at 0.4 m s−1 for a week before the onset period and, during that time, the other ADSM parameters and M did not fluctuate significantly. At the beginning of the onset period, glacial transit speed increased to 0.65 m s−1, and the fraction of mobile water, β, increased, whereas the fraction of returned tracer mass, M, dropped from 0.8 to 0.4. Towards the end of the onset period, lake discharge had increased to 3 m3 s−1 and glacial transit speed to 0.95 m s−1, accompanied by a sharp drop in dispersion, D. The lake discharge reached ∼5.5 m3 s−1 the day after the onset period. At the same time the glacial transit speed reached its maximum of 1.05 m s−1 and then decreased to values in the range 0.9–0.7 m s−1, typical for a well-developed moulin on Gornergletscher (Reference WerderWerder, 2009).
Figure 3 shows the subglacial water pressure heads of five boreholes. BH6, located 50 m from M4, responded strongly to the lake drainage: its water level stayed high throughout the onset period, and afterwards it fluctuated, again diurnally. BH4, 400 m south and upstream of M4, also responded strongly. Note that the decrease in water pressure was delayed by 8 hours, compared to BH6. The water pressure in boreholes BH1, BH2 and BH3 did not react to the lake drainage. In contrast, when the lake drained as a subglacial outburst flood, a response was seen for several days in all these boreholes (Reference Huss, Bauder, Werder, Funk and HockHuss and others, 2007; Reference Werder, Loye and FunkWerder and others, 2009).
Fig. 3. Subglacial water pressure head, h, measured in boreholes in the vicinity of the lake and connection channel. The two vertical lines delimit the onset period.
Figure 4 presents complementary measurements around the onset period. The water exiting the lake was 2°C on the first day of the onset period, then dropped to 0.8°C during the night and rose to 1°C on the second day of the onset period (Fig. 4b). Figure 4d shows subglacial water pressure head data from BH1 in 2006 (black curve, also shown in Fig. 3e) and data, also from BH1 (grey curve), for the year 2005. The peak-to-peak amplitude of the diurnal water pressure head fluctuations in BH1 were 25 m in 2006 and 150 m in 2005. Note that the water pressure head in 2006 stayed level at ∼325 m during the time of low pressure and then increased during the day. The rain events on 5 and 6 July (Fig. 4c) can be recognized in the proglacial discharge and in BH1 and BH2. The colder temperature on 7 July led to lower proglacial discharge and a drop in water pressure in BH1 and BH3.
Fig. 4. (a) Proglacial discharge (Q, black, left scale) and lake discharge (grey, right scale); (b) temperature of water entering M4, T water; (c) air temperature (T air, black) and precipitation (P, grey); and (d) pressure head in borehole BH1 (black curve shows data from 2006, grey curve data from 2005). The two vertical lines delimit the onset period.
Jökulhlaup Model
We compare flow speeds inferred from tracer experiments with water flow speeds calculated with Reference ClarkeClarke’s (2003) model. So far this model has been calibrated using only lake and proglacial discharge hydrographs. Thus this study presents the first more stringent test of the model, comparing it to inferred water flow speeds and not only to the lake discharge.
Model formulation
Reference ClarkeClarke’s (2003) jökulhlaup model integrates a modified version of the Spring–Hutter equations (Reference Spring and HutterSpring and Hutter, 1982) which describe water flow through R channels. These R channels (Reference RöthlisbergerRöthlisberger, 1972) are en- or subglacial circular or semicircular channels incised into the ice. These equations describe the time evolution of water pressure, channel cross-section, water flow speed and water temperature. The Spring–Hutter equations contain several empirical relations, two of which are of interest here: firstly, the relation between pressure gradient and water flow speed and, secondly, the heat transfer relation.
The Manning–Gauckler–Strickler formula relates water flow speed, v, to hydraulic head gradient, dH/ds:
where n man is the Manning roughness, s is the along-channel coordinate and R h is the hydraulic radius defined as
where S is the cross-sectional area of the channel and P w is the wetted perimeter. The Manning roughness is one of the tuning parameters of the model. The same average roughness is used for englacial and subglacial channels. Reasonable values of n man lie in the range 0.02 < n man < 0.08 m−1/3 s (e.g. Reference ClarkeClarke, 2003).
The amount of melted ice, m, at the conduit walls per unit length is proportional to the temperature difference between ice and water, ΔT:
where L is the latent heat of fusion and U h is the heat transfer coefficient. The model uses the following empirical relation (Reference McAdamsMcAdams, 1951):
where Re is the Reynolds number, Pr is the Prandtl number, κ is the thermal conductivity of water and P m is the ice-walled perimeter of the conduit. (Reference Bird, Stewart and LightfootBird and others, 1960, give a list of alternative heat transfer relations.)
As input the model needs the geometry of the connection channel and the thickness of the overlying ice. Channels of different sinuosities, σ, are implemented in the model by stretching the along-channel coordinate, s, by σ. The boundary conditions at the inlet are the water level and the water temperature. At the outlet the water pressure is needed as a boundary condition. Reference ClarkeClarke’s (2003) code has been modified such that the water level at the outlet (normally used to simulate a jökulhlaup which terminates underwater) can be a function of time. The fitting parameters are the initial cross-sectional area of the connection channel, c 0, its sinuosity, σ, and its roughness, n man . The model was fitted to the lake discharge and for some runs also to the inferred water flow speeds using, first, a grid search on the tuning parameters (n man , σ, c 0) and then a steepest-descent method.
Connection channel location
Only the connection channel, connecting the lake to the main drainage channel of the glacier, is subjected to a jökulhlaup-like evolution (grey arrow in Fig. 1) and hence can be simulated using Reference ClarkeClarke’s (2003) model. We consider a connection channel transit path lying in the region of the grey shaded area in Figure 1. We restrict it to this area since tracer experiments from both M2 and M1 (Fig. 1) show only a modest response to the subglacial outburst floods in 2005 and 2007 (Reference Werder, Loye and FunkWerder and others, 2009), hence connection channels lying north of M2 are not considered. The response of BH4 during the onset period can be explained by a flow direction in the lake dam area parallel to the crevasse orientation, as was also determined by tracer experiments (Reference WerderWerder, 2009). The direction of the crevasses is north–south and BH4 lies upstream of M4. BH3 did not show a response either, so connection channels south of BH3 are not considered, which would also grossly disagree with the direction of steepest descent of the hydraulic potential. Borehole BH2 also does not show a response, which is puzzling as it lies in the middle of the region of possible transit paths. Active source seismology (conducted by lowering explosives into M4) showed that the moulin/channel shaft drops very steeply (crosses on Figs 1 and 5 mark shot locations) and is orientated slightly southwards, suggesting that the connection channel lies south of BH2. A connection channel lying in the grey shaded area (Fig. 1) has a transit distance of 
.
Fig. 5. A cross-section of the glacier with the assumed connection channel in grey. The crosses indicate locations of active seismic shots conducted inside M4. The boundary conditions needed at the inlet and outlet are given in parentheses.
The geometry of the channel in the vertical is depicted in Figure 5, which shows a cross-section of the glacier along the connection channel. We assume that the channel drops steeply englacially, as supported by the active source seismology, to reach the bed and then follows it. The length along the connection channel, having a transit distance 
, is l
con
= 790 m. The vertical geometry of channels lying in the grey shaded area is very similar, and thus their geometry can be adjusted to different transit distances by stretching it in the horizontal. This leads to a range of 790 < l
con
< 1120 m. To finally arrive at the channel flow path length, we assume the connection channel sinuosity to be constrained by 1 ≤ σ < 2. We absorb the above range of l
con
into the sinuosity, fix l
con
= 790 m and vary 1 < σ < 2.8, of which a factor of up to 1.4 can be due to different connection channel transit paths.
Flow speed estimate
An estimate must be made of the tracer flow speed in the connection channel in order to be able to compare the water flow speeds from the model with the experimental data. We divide the total residence time, 
, i.e. the time the tracer took in experiment i to reach the detection station, into the sum of three parts: the proglacial residence time, 
, the main channel residence time, 
, and the connection channel residence time, 
:
where 
was measured and 
can be determined from Equation (1). 
is unknown but must be greater than 0 to assure a finite flow speed in the main channel. To get an upper bound on 
, we assume that 
does not change much due to the additional lake water influx or due to changing meltwater discharge conditions from day to day, as suggested by the constant amplitude and mean of the proglacial discharge (Fig. 4a). Hence an upper bound for 
is the shortest measured glacial residence time
where experiments conducted in the morning, early afternoon and evening are treated separately, i.e. a different t shortest is obtained for each time of day. All of them achieved the shortest glacial residence time on 7 July, the day after the onset period. The t shortest for morning, early afternoon and evening correspond to main channel transit speeds of 0.65, 0.86 and 0.79 m s− 1, respectively.
The inequality on 
,
together with Equation (6) leads to upper and lower bounds for 
:
These bounds on 
rest only on the assumption that 
does not exceed t
shortest. To translate this into flow speed bounds which can be compared to Reference ClarkeClarke’s (2003) model, an assumption on the connection channel flow path length needs to be made, as discussed in the previous subsection:
The values of 
will be referred to as inferred flow speeds or flow speed bounds. Note that the flow speed bounds depend on the sinuosity, σ.
Model set-up
The model is run with the channel geometry discussed above (Fig. 5). The channel cross-section is set to circular for the englacial part and to semicircular for the subglacial part. To reduce the number of free parameters, we use an average Manning roughness and do not distinguish between the roughness of the ice walls and the glacier bed. The boundary conditions at the channel inlet are the measured water temperature (Fig. 4b) and lake level, which was constant during the onset period. At the channel outlet we use the water pressure measured in borehole BH1 and also data from the same borehole but from the previous year, 2005, as the borehole was not so well connected in 2006 (Fig. 4d).
Four different model runs (Mod1–Mod4) are used to investigate the influence of different boundary conditions, physics and fitting procedures. We tune model Mod1 to the lake discharge only and Mod2–Mod4 are tuned, in addition, to the inferred flow speeds. The tuning parameters are the Manning roughness, n man , the initial channel size, c 0, and for Mod2–Mod4 also the sinuosity, σ.
It is found that the modelled flow speeds in Mod1 are higher than the inferred flow speeds. If the sinuosity is increased, the inferred flow speed also increases, which is simulated in Mod2 and Mod3. Another possibility is to lower the modelled flow speeds. For pressurized channel flow, the discharge is equal to the flow speed times the channel cross-sectional area. Hence, to decrease the flow speed for a given discharge, the cross-sectional area needs to be larger. This can be achieved by increasing melt of the channel walls (i.e. the heat transfer) or by reducing creep closure. Creep closure can be neglected in these settings as the high water pressure prevents it. This leaves the heat transfer, which is increased by a factor of 2 in Mod4. Table 1 gives a summary of the four model runs performed.
Table 1. Summary of the settings and results (n man , σ and c 0) for model runs Mod1–Mod4
Model Results
Mod1 (Fig. 6a and b) is fitted to the lake discharge using n man and c 0 as tuning parameters and setting σ = 1. The roughness is high at n man = 0.071 m−1/3 s. The lake discharge fits the measurements well, but the model speeds are ∼0.55 and 0.25 m s−1 too high on the first and second day of the onset period, respectively.
Fig. 6. Results of model runs Mod1–Mod4. Top panels show the modelled discharge, Q (solid curve) and measured discharge (crosses). Bottom panels show modelled flow speeds, v con (solid curves) and inferred flow speed bounds (bars). The legend gives the name, the Manning roughness, n man , the sinuosity, σ, and the initial channel cross-section, c 0, of each model run. The dashed lines delimit the onset period. (a, b) Mod1 fitting only the discharge; (c, d) Mod2 also fitting the inferred flow speeds; (e, f) Mod3 using subglacial water pressure data from 2005 with higher diurnal variations; (g, h) Mod4 using enhanced heat transfer (see Table 1).
Mod2 (Fig. 6c and d) is fitted to both Q and v con using all three tuning parameters (n man , c 0 and σ). The model fits the inferred flow speeds of the second day well, but not those of the first day (0.5 m s−1 difference). Mod2 has, compared to Mod1, a larger sinuosity and roughness of 1.43 and 0.049 m−1/3 s, respectively. Note that, due to the larger sinuosity, the upper flow-speed bounds estimated by Equation (10) become greater than for Mod1.
Mod3 (Fig. 6e and f) is fitted as Mod2 but uses 2005 subglacial water pressure data as a lower boundary condition. The modelled discharge has larger diurnal fluctuations and fits less well; in particular, its increase at the end of the onset period is far steeper than that of the measurements. Also, the pronounced modelled flow speed peak in the morning of the second day cannot be seen in the inferred flow speeds. The modelled flow speeds are 0.4 m s−1 too high on the first day and lie within the flow-speed bounds on the second day. The roughness is comparable to Mod2, whereas the sinuosity is greater (1.86).
Mod4 (Fig. 6g and h) is run with a two-fold increased heat transfer. It is fitted to both Q and v con using n man , σ and c 0. The modelled flow speeds are lower than in the other models and fit well on the second day but are still 0.3 m s−1 too high on the first day. The Manning roughness, n man = 0.105 m−1/3 s, is large and the sinuosity is fitted to 1.00.
The time evolution to 1.5 days beyond the onset period of the modelled discharges, flow speeds and outlet water temperatures is shown in Figure 7. The modelled discharges, Q, diverge and eventually reach 20–28 m3 s−1. The modelled flow speeds, v con , are between 2.2 (Mod4) and 3.4 m s−1 (Mod2). The connection channel outlet water temperatures are all within 0.1°C of each other. The output of Mod3 exhibits much larger fluctuations, due to the larger pressure fluctuation at its outlet. Note that these model results cannot be compared to measurements, as the lake drainage was spillway-limited after the onset period. This means that the water flow inside the moulin was partly open-channel flow whereas the model only works for fully pressurized flow conditions.
Fig. 7. The simulated jökulhlaup evolution up to 1.5 days beyond the onset period of Mod1–Mod4. (a) Lake discharge, Q (measured discharges are marked with crosses); (b) flow speed, v con ; (c) water temperature at the connection channel outlet, T.
Comparison of measured lake discharge data and inferred tracer flow speeds to results obtained with Reference ClarkeClarke’s (2003) model shows that:
Fitting the model only to the lake discharge using the shortest possible connection channel yields model speeds that are too large and a high Manning roughness (Mod1).
The model can be fitted to the inferred flow speeds on the second day with a higher sinuosity (or a longer connection channel).
Alternatively, the model speeds can be fitted to the inferred flow speeds by enhancing the heat transfer. This leads to a very high Manning roughness (Mod4).
It is not possible to discriminate between Mod2 and Mod4 on the basis of the presented experimental data.
The low inferred flow speeds on the first day cannot be modelled by any parameter combination that would give a reasonable fit on the second day.
All the modelled flow speeds are either higher than or near the upper end of the experimental bounds.
Mod3, using the subglacial water pressure recorded in 2005 as the lower boundary condition, results in a poorer fit than Mod2 or Mod4.
Discussion
Observations
The tracer experiments show a pronounced evolution of the drainage path of M4 during the onset period (Fig. 2). The increase in glacial transit speed on the first day of the onset period is mainly due to the increased pressure difference driving the water. The further increase on the second day is caused by the enlargement of the channel. Before the lake drained, the channel was inefficient but then developed quickly, as can be seen from the reduction of dispersion, D, and the increased fraction of mobile water, β. However, contradictory to this are the large diurnal water-level fluctuations in BH6 (located only 50 m from M4), which are commonly assumed to be caused by an efficient drainage system. The trace right at the beginning of the onset period shows a greatly reduced fraction of returned tracer mass, M. This might be caused by water being pushed out laterally from the channel into storage areas. On the second day of the onset period the local storage areas may have been filled and the fraction of returned tracer mass increases again to previous values. The data from BH4 (Fig. 3b) also indicate the release of stored water: its water pressure drops 8 hours after the pressure in BH6 drops, which could be due to water flowing out of storage keeping the pressure high. This 8 hour time lag in the signal of BH4 measurements is also present before and after the drainage, indicating that storage-release processes also play a role during normal discharge regimes (cf. Reference Hock, Iken and WanglerHock and others, 1999). The time lag also indicates that BH4 lies upstream of BH6, because otherwise the water pressure in BH6 would lag the pressure in BH4. Hence, water flow does not always follow the direction of the gradient of the calculated hydraulic potential (Fig. 1; also Reference WerderWerder, 2009).
The bounds on the flow speed in the connection channel (Fig. 6) depend on the sinuosity, σ, of the channel and range from a maximum speed of 1.4 m s−1 for σ = 1 to 2.7 m s−1 for σ = 1.86. The inferred speeds on the first day of the onset period are slow (<0.6 m s−1). Combined with the high dispersion, this suggests that the drainage system was not very efficient, as could be the case for a braided flow path or flow through a network of englacial cracks (Reference Fountain, Jacobel, Schlichting and JanssonFountain and others, 2005). The lower bounds on the flow speeds are always very low due to the poor constraint on the upper bound of 
, but from the comparison to the model, it seems that the upper bounds on the flow speed are more relevant.
Comparison measurements and model
Even though the lake drained by overspilling into a moulin on its shore, the situation was similar to a subglacial outburst flood during the onset period. The main difference is that during a subglacial outburst, the water pressure is higher at the channel inlet which is then located at the lake bottom. However, this difference makes the used model (Reference ClarkeClarke, 2003) no less applicable, as the upper boundary condition can be adjusted accordingly and thus we can test it against our measurements.
Mod1 (Fig. 6a, b) shows that a good fit to the lake discharge can be achieved without matching the flow speeds inferred from the tracer experiments. Hence, discharge measurements alone are not enough to test the Spring–Hutter equations applied to a jökulhlaup. In order to match the inferred flow speeds, either the sinuosity/length of the connection channel (Mod2, Mod3) or the heat transfer (Mod4) need to be increased. However, on the first day of the onset period, the model cannot fit the inferred flow speeds while still matching the lake discharge and the inferred flow speeds of the second day. Thus, there is a qualitative difference in the drainage path between the first and the second day of the onset period. This is supported by the other observations discussed above, in particular the decrease in tracer dispersion. The inability to model the flow speeds on the first day of the onset period may be due to the following reasons: the water does not initially flow through R channels (e.g. braided channels or englacial cracks); the Spring–Hutter equations are not valid for R channels at these low flow speeds and discharges; the geometry changes profoundly between the first and second day (e.g. the sinuosity decreases considerably). Any of these reasons indicates that an accurate model of the very beginning of a jökulhlaup needs to go beyond the Spring–Hutter equations. So far only one jökulhlaup model does this and includes also sheet flow (Reference Flowers, Björnsson, Pálsson and ClarkeFlowers and others, 2004), which probably could capture flow through a more distributed system with suitably chosen parameters.
On the second day of the onset period the model performs well if either the sinuosity/channel length is increased or the heat transfer from the water to the channel wall is increased. A more sinuous or longer channel (Mod2, Mod3) leads to a reduced roughness and higher flow speeds to fit the discharge. However, because the inferred flow speed also depends on the sinuosity of the channel, the model results match up with the measurements. For Mod4, the heat transfer was increased, as it has been noted that the outlet water temperature is overestimated by the Spring–Hutter equations (Reference JóhannessonJóhannesson, 2002; Reference ClarkeClarke, 2003); however, Figure 7c shows that this does not significantly reduce the water temperature. The use of a two-fold enhanced heat transfer coefficient (Equation (4)) leads to a very high Manning roughness (n
man
= 0.105 m−1/3 s) to keep flow speeds low in the larger channel. A high resistance to water flow can be caused by a channel which is low and broad. Reference Hooke, Laumann and KohlerHooke and others (1990) introduced low and broad channels with the shape of a circular segment characterized by its central angle, θ. Such a Hooke channel has a circular cross-section for θ = 360° and a semicircular one for θ = 180°. The hydraulic radius, R
h (Equation (3)), of a Hooke channel becomes smaller with decreasing θ for a fixed cross-sectional area, S. The hydraulic radius in terms of S and θ is given by 
. A smaller hydraulic radius reduces water flow speed (Equation (2)). Hence, the high Manning roughness of Mod4 can be reduced to an effective value of n
man
= 0.049 m−1/3 s (like Mod2) by the smaller hydraulic radius of a Hooke channel with a central angle of θ = 16.5° (this geometry also leads to higher closure rates, but these are not important in these settings). This channel shape is very low and broad but could occur during the onset of a jökulhlaup if the drainage path evolved from a braided one, as suggested above.
Models Mod1, Mod2 and Mod4 produce virtually identical discharge hydrographs and are distinguished only by their water flow speeds. Hence, the comparison to measured flow speeds is needed for further model validation and leaves Mod2 and Mod4 as options. Unfortunately, it is impossible to discriminate between these two using our experimental data. This is because the flow speeds inferred from the tracer experiments are themselves dependent on the sinuosity. The water temperature at the outlet is almost identical (Fig. 7c) for Mod2 and Mod4, indicating that it is difficult to distinguish high sinuosity from enhanced heat transfer solely by temperature measurements. However, in situations where there are no direct flow speed measurements, low water temperatures at the outlet can indicate high sinuosity or enhanced heat transfer since Mod1 produces higher outlet temperatures than either Mod2 or Mod4. A model of the whole subglacial drainage system would not need the borehole water pressure data as a boundary condition, thus allowing a comparison for validation. However, many more assumptions would be needed to set up such a model and thus it is unlikely to yield more accurate results.
The presented interpretations are based on a model which was deliberately kept simple; nonetheless, there are a few potential shortcomings. An argument was presented to constrain the location of the connection channel to the grey area in Figure 1. However, it is also plausible, from the hydraulic potential, that the connection channel could pass between M1 and M2. Also, the connection channel is likely to have tributaries, influencing its discharge, initial size and pressure. This would lead to increased tracer flow speeds, because in R channels the flow speed is an increasing function of discharge. Hence, the model, not taking into account tributaries, should underestimate rather than overestimate water flow speeds. The water pressure boundary condition at the connection channel outlet is insufficiently known. Figure 6e and f show that the use of the previous year’s subglacial water pressure data with large diurnal fluctuations as a boundary condition leads to variations in the modelled speeds and discharge that are greater than those observed in the measurements. In all likelihood, the diurnal fluctuations found in the pressure data recorded during the onset period were too small, whereas those from the previous year’s data were too great (Fig. 4d). The additional water coming from the lake would keep the subglacial water pressure higher during the night, as was observed, for example, during the 2005 outburst on 11 and 12 June (Reference Werder, Loye and FunkWerder and others, 2009), when the nightly pressure minima became less pronounced at lake discharges comparable to those encountered during the presented measurements. This higher subglacial water pressure at night would suppress the large modelled variations. This sensitivity on the pressure boundary condition at the connection channel outlet shows that the rest of the subglacial drainage system, dictating the subglacial water pressure, can have a great influence on the trigger and initial jökulhlaup evolution. Hence, to simulate the onset of a lake drainage, the prevailing conditions in the subglacial drainage system need to be taken into account. The model does not include storage–release processes, for which there is clear evidence from tracer experiments and borehole water pressure data (see also Reference Huss, Bauder, Werder, Funk and HockHuss and others, 2007). They can influence the pressure and discharge conditions considerably, in particular during the onset. Water diverted into storage should lower the tracer flow speed by the fraction of diverted water. The inferred flow speeds obtained from the tracer measurements on the first day are about two to three times less than the modelled speeds. This could be caused by forcing one-half to two-thirds of the water into storage during this time. On the first day, tracer return mass diminished by half, suggesting that the low flow speed could be caused for the most part by storage processes.
Conclusions
For the first time, water flow speeds calculated by a jökulhlaup model were compared to speeds inferred from dye-tracer injections using the lake outlet. This showed that an ingenuous fitting to the discharge led to model flow speeds that were too high. During the first day of the lake drainage the inferred flow speeds were too low to be fitted at all, indicating that the water was not yet flowing in an R channel. Hence, to simulate the onset of a jökulhlaup it is necessary to go beyond the Spring–Hutter equations. Afterwards, the inferred flow speeds can be fitted by either increasing the sinuosity of the connection channel or by enhancing the heat transfer. For enhanced heat transfer, the roughness needs to be increased to an improbably high value for a semicircular channel, but it can be reduced to reasonable values assuming a low and broad channel. However, the experimental data do not allow us to discriminate between these two possibilities. In particular, a direct measurement of the location and sinuosity of en-and subglacial channels would be desirable, as would a better understanding of the heat transfer in an R channel. The former seems difficult to achieve, but the latter could be achieved by laboratory experiments.
Acknowledgements
The project was funded by the Swiss National Science Foundation, grant Nos. 200021-103882/1 and 200020-111892/1. We are grateful to many members of VAW/ETH Zürich and others who helped with data acquisition, fieldwork and data processing. Grande Dixence SA allowed us use their infrastructure to install the fluorometers and provided the discharge data. The Swiss military and the International Foundation High Altitude Research Stations Jungfraujoch and Gornergrat (HFSJG) provided logistical support. M. Schurter lent us his BackScat fluorometer. G.K.C. Clarke provided his computer program to integrate the Spring–Hutter equations. M. Huss calculated the lake discharge. H. Horgan, P. Riesen, M. Truffer, M. Huss and, in particular, T.V. Schuler made many helpful comments and suggestions. We also thank M.J. Roberts, G.K.C. Clarke and the scientific editor R. Hock for their numerous and helpful revisions.
