Glaciological data

Surface mass-balance measurements (Reference Østrem and BrugmanØstrem and Brugman, 1991) are acquired using 6 m stakes drilled into the ice/firn where the length of the exposed stake is measured on an annual/seasonal basis. Additional accumulation measurements are acquired across the glacier by snow depth sounding or radar profiling, and the snow density is measured at selected stakes in April/May. In this study, only measurements at stake locations are used for model calibration because they are consistent through the entire time series. Summer measurements are typically acquired in early September. Due to weather-related inaccessibility, not all stakes were measured at the end of summer every year. Therefore, we use the net change in exposed stake heights between two successive April/May measurements minus the previous winter’s accumulation to derive the summer balance. In situ mass-balance data are available from 1987 to 2008 for Kongsvegen and from 2003 to 2010 for Kronebreen. Measurement details are given by Reference Melvold and HagenMelvold and Hagen (1998), Reference Hagen, Melvold, Eiken, Isaksson and LefauconnierHagen and others (1999) and Reference BaumbergerBaumberger (2008). Due to the heavily crevassed glacier tongue of Kronebreen, mass-balance measurements were not obtained below 500ma.s.l. before 2008. Three stakes in this region have been maintained since 2008. However, potential winter ablation has not been successfully measured, due to the disappearance of the stakes between autumn and spring.

Meteorological data

The Norwegian Meteorological Institute has measured precipitation and temperature in Ny-Ålesund since 1969. In 1974 the weather station was moved to ~34 m below its original location. Small lapse corrections (0.2°C) have been shown to have little significance in correlations with mass balance (Reference Lefauconnier and HagenLefauconnier and Hagen, 1990). Synoptic weather observations were carried out (three or four times daily), until 1994 when an automatic weather station (AWS) began continuously monitoring at an hourly recording interval (Norwegian Meteorological Institute; http://www.met.no). Before 1994, precipitation measurements were collected manually using a rain gauge; since 1994 the weight of the collected precipitation in the rain gauge has been measured hourly (http://www.met.no). We use daily accumulated precipitation and daily average and maximum temperatures to drive the surface mass-balance model. Average and maximum temperatures are employed to determine whether partial-melting days (i.e. days in which melt occurs but the average temperature is below freezing) influence the model parameter values. In addition, daily average data from AWSs on Kongsvegen and Kronebreen are used to validate distribution of temperature data from Ny-Ålesund (Fig. 1).

Elevation data

The elevation data consist of DEMs and a GPS profile (Table 1). The 1966, 1990 and 1995 DEMs were generated from vertical aerial imagery acquired in summer, made on a digital photogrammetric workstation (Reference AltenaAltena, 2008). The 2007 DeM is a SPOT5-HRS (high-resolution stereo) product (GES07-043-NorthWestSvalbard) obtained through the IPY- SPIRIT project (International Polar Year SPOT5 stereoscopic survey of Polar Ice: Reference Images and Topographies; Reference Korona, Berthier, Bernard, Rrsmy and ThouvenotKorona and others, 2009). The DEM is created from satellite stereo imagery of a forward- and backward-looking sensor with a base-to-height ratio of 0.8 (Reference Bouillon, Bernard, Gigord, Orsoni, Rudowski and BaudoinBouillon and others, 2006). In Svalbard, SPOT5-HRS DEMs are reported to have an elevation accuracy of 3 m and a precision of 5 m relative to the Ice, Cloud and land Elevation Satellite (ICESat) (Reference Korona, Berthier, Bernard, Rrsmy and ThouvenotKorona and others, 2009; Reference Nuth and KääbNuth and Kääb, 2011).

Since the DEMs from the 1990s do not cover both glaciers fully in a single acquisition, we use a 1995 DEM for Kongsvegen and a 1990 DEM for Kronebreen. The data in 1966 and 1990 were incomplete in some regions, due to low visible contrast (i.e. higher elevations), rendering image matching difficult. The 1966 DEM gaps are filled with the 1966 contour data made previously on an analogue photogrammetric workstation using the same images. The 1990 data on Kronebreen are missing above 700ma.s.l. (40% of the glacier area) which we fill using a 1996 differential GNSS elevation profile measured along the approximate center line every 10s (~40m). The 1996 profile data are elevationally adjusted to 1990 by removing the mean elevation difference in an area of overlap and assuming the glacier slope has not changed. The 1995 DEM of Kongsvegen is not missing any data.

Mass continuity

Glacier-elevation changes, mass balance and flux of ice in and out of the system are related through the equation of mass continuity. A detailed derivation and discussion of glacier mass continuity is given by Reference Cuffey and PatersonCuffey and Paterson (2010). Here we simplify the discussion and list the assumptions related to combining elevation changes and the surface mass balance. For any volume element of a given material having density ρ, mass continuity is written:

where is a velocity vector and β is a production term. Eqn (1) states that any local change in density is balanced by the net flux of material into or out of the considered volume plus any source or sink of mass. Integration from the bed, h _b, to the surface, h _s, provides the continuity equation for a vertical column through a glacier:

Where

The term on the left-hand side of Eqn (2) represents the change in mass of the given volume and may be approximated through gravity variations (e.g. Reference Wahr, Swensen, Zlotnicki and VelicognaWahr and others, 2004; Reference Luthcke, Arendt, Rowlands, McCarthy and LarsenLuthcke and others, 2008) or elevation changes. On the right-hand side, the mass balance, b, is the sum of surface (b _s), englacial (b _e) and basal (b _b) components (Eqn (3)). Eqn (4) defines the horizontal flux, , as the column average velocity, , multiplied by the density, ρ, which closely approximates the density of ice in cases where firn is a small proportion of the ice thickness (Reference Cuffey and PatersonCuffey and Paterson, 2010).

To relate the mass change (left-hand side of Eqn (2)) to observed elevation changes, ∂h/∂t, we introduce an effective density, ρ _eff, to express the temporal and vertical changes of the column density:

We continue in water-equivalent units because this is the common measuring practice for glacier mass balance, b. We now assume that englacial and basal mass balances are negligible (b _e ≪ b _s and b _b ≪ b _s) and that the flux divergence, , represents mass change of incompressible ice. The continuity expression can then be reduced to

where ρ_w is the density of water and k is a conversion factor from height differences to water-equivalent changes:

If ∂h/∂t is observed over a significantly long time period, k can be approximated by the density ratio of ice to water (0.9) below the equilibrium-line altitude (ELA). This is because small changes in the less dense snow have little impact on the column average density and thus all changes will be of incompressible glacier ice. In the firn area, changes in the proportion of firn in the ice column can alter k, due to the compressibility of firn. It is often assumed that firn thickness and density are constant through time (Sorge’s law; Reference BaderBader, 1954), in which case k = 0.9.

Solving mass continuity over the entire glacier system requires integration of Eqn (6) over the glacier surface area, A:

with all terms in water equivalent. Eqn (8) relates the volume change, ∂V/∂t, to the water-equivalent flux divergence, , and the glacier-wide mass balance, B. Applying the divergence theorem to the last term in Eqn (8) results in the relationship between the glacier-wide integrated flux and the water-equivalent flux through a boundary, R

and substitution into Eqn (8) results in

where is the normal vector to the closed boundary, R. The second term on the right-hand side may represent the influx of ice by avalanching or the loss through calving and potential submarine melting on the calving face.

Using Eqn (10), we consider two solutions: noncalving and calving glaciers. Often for non-calving glaciers, is assumed equal to zero, resulting in

This has formed the basis of many comparison studies aimed at detecting systematic errors in the cumulative, glacier-wide in situ surface mass-balance values by using geodetically measured volume changes (e.g. Reference KrimmelKrimmel, 1999; Reference Elsberg, Harrison, Echelmeyer and KrimmelElsberg and others, 2001; Reference Cox and MarchCox and March, 2004; Reference Thibert, Blanc, Vincent and EckertThibert and others, 2008; Reference ZempZemp and others, 2010). For a calving glacier, is equal to the flux, Q, through the boundary, R, of the glacier:

Practically, solving mass continuity of an entire glacier system requires definition of the boundary geometry. For simplicity, or due to lack of updated maps, this geometry may be held constant. For example, Reference Elsberg, Harrison, Echelmeyer and KrimmelElsberg and others (2001) introduce the concepts of a ‘reference’ and ‘conventional’ surface for mass-balance integration and suggest a transformation between them. The ‘reference’ surface is a constant map (elevation and extent) and is considered to be more climate-related, as it ignores the effects of surface change on the mass balance. The ‘conventional’ mass- balance surface is the actual mass change of the glacier relevant for hydrological and sea-level change studies. In the following, we consider mainly the change in glacier extent through time due to the stronger influence on mean mass balance than that of a changing surface elevation (Reference PaulPaul, 2010). Nevertheless, changing surface elevations are roughly compensated for in the surface mass-balance model through our use of an average elevation of the two DEMs within an epoch. Eqn (12) can be modified to handle the volume of retreat/advance separately:

where

Derived in this way, ∂V/∂t and B of Eqns (11) and (12) can be solved using a ‘reference’ surface, defined as the smallest glacier area. Q is then the ice flux through the crosssectional area (flux gate) defined by the glacier front at the time of smallest glacier area (i.e. the most recent area in cases of retreat). In Eqn (13), ∂V _r/a/∂t, and B _r/a are the volume change and mass balance of the receding (r) or advancing (a) area. Quantifying ∂V _r/a/∂t is not a problem, provided we have knowledge of the basal elevation or depth below sea level in the retreat or advance area. B _r/a can be solved by assuming a linear retreat of the front. Q _r/a is defined as the retreat/advance flux and is the difference between the net flux into or out of the glacier, Q′, and the flux out of the gate defined by the front of the smallest area, Q. Hence, Q _a > 0 and Q _r < 0.

Elevation and volume changes

Elevation changes are calculated for two time epochs using the data listed in Table 1: 1966-1995-2007 for Kongsvegen and 1966-1990-2007 for Kronebreen. The 1966, 1990 and 1995 DEMs are resampled to 40 m (SPOT5-HRS DEM resolution) using a 4 × 4 window average block filter on the original 10m products. DEM pairs are co-registered using the sinusoidal dependency between the vertical differences over stable terrain with aspect (Reference Nuth and KääbNuth and Kääb, 2011). The most distinct shift between the products occurs in the 2007 SPOT5-HRS DEM because the 1966 and 1990 data are obtained using roughly the same ground-control points. The triangulation of shift vectors between three DEMs for Kongsvegen and Kronebreen shows a horizontal and vertical coherency of <2 m (Table 2). This vertical residual represents the accuracy after removal of vertical systematic biases, and becomes noticeable when analyzing the mean glacier elevation change. The stochastic errors (precision) of the datasets are generally better than 10 m (σ in Table 2). No elevation-dependent biases are detected.

Estimating volume changes from elevation changes can be accomplished using two methods, typically depending on whether the data available are center-line profiles (Reference Arendt, Echelmeyer, Harrison, Lingle and ValentineArendt and others, 2002) or full-range DEMs (Reference Berthier, Schiefer, Clarke, Menounos and RemyBerthier and others, 2010). The volume change derived from differential DEMs is the summation of the elevation-change pixels over the glacier surface multiplied by the pixel area (‘grid’ method). If only center-line data are available, a ‘hypsometric’ approach can be employed, assuming elevation changes across an elevation bin are normally distributed (Reference Berthier, Arnaud, Baratoux, Vincent and RemyBerthier and others, 2004). In this method, volume change is derived by summation of the product between elevation change as a function of elevation (average or centerline assumed average) and the glacier hypsometry (area- altitude distribution). Converting volume change to mass change requires knowledge of the effective density, typically assuming all elevation changes consist of changes in ice thickness ( k = 0.9).

Surface mass-balance (SMB) model

The SMB model is driven at a daily time-step using meteorological data from Ny-Ålesund over the period 19692007. The model is calibrated for each glacier individually using the entire time series of available point mass balances (Table 1) in order to derive model parameters individually for each glacier. This decision is based upon a priori knowledge that Kongsvegen receives, on average, slightly more accumulation at its highest elevation (800 m) than Kronebreen at 1400 m and due to the large superimposed ice zone on Kongsvegen (Reference Kbnig, Wadham, Winther, Kohler and NuttallKctnig and others, 2002; Reference LangleyLangley and others, 2007; Reference Brandt, Kohler and LuthjeBrandt and others, 2008) and the possible lack of one on Kronebreen (Reference BaumbergerBaumberger, 2008).

Geometry changes

By 2007, Kongsvegen and Kronebreen had lost about 5% and 1% of their 1966 area, respectively. The majority of area loss occurred in the 1966-1990/95 epoch and the front position has remained relatively stable since 1990/95. Within the time frame of this study, the medial moraine dividing the two glaciers has migrated, reflecting the mass flux variations between the two glaciers (Reference Kääb, Lefauconnier and MelvoldKääb and others, 2005). In 1966 Kongsvegen occupied ~50% of the calving front, whereas in 2007 the glacier rested on land. The elevation changes for the different epochs are shown in Figure 5. The retreat of the glacier front is apparent (Fig. 5f and g) as increasingly negative elevation changes between 0 and 150 m a.s.l. Retreat from 1990 to 2007 is negligible.

Fig. 5. Elevation-change rates on Kongsvegen, Kronebreen and Kongsbreen for (a) the entire time series, 1966-2007, (b, d) epoch I (1966— 1990/95) and (c, e) epoch II (1990/95—2007). Elevation losses over the marine retreat area do not include ice below sea level. The gray solid lines in (a) represent the center lines of both glaciers. (f, g) Elevation-change rates as a function of elevation for each glacier. In gray are all ∂h/∂t pixels from the entire time series (1966—2007) where full spatial distribution is possible. ∂h/∂t for epochs I and II are shown for the center line (solid circles), GPS (unfilled circles) and 50 m elevation interval averages (unfilled squares).

On Kongsvegen, the elevation change rate remained relatively stable between the two time epochs (Fig. 5f). Above 550 m a.s.l., the glacier surface elevation increases, greatest along the center line (up to +1 ma^-1) compared with the average of the elevation interval. This is especially true for epoch I, less so during epoch II, which is in agreement with the decreased thickening observed by Reference Hagen, Eiken, Kohler and MelvoldHagen and others (2005). Also, elevation increases occur lower (~400 m) on Sidevegen than on Kongsvegen (~500m), visible as the cloud of gray dots between 300 and 500 m a.s.l. in Figure 5f that deviate from the Kongsvegen center line.

On Kronebreen, ∂h/∂t of the glacier tongue for the full period (1966-2007) is similar to Kongsvegen as well as Kongsbreen (Fig. 1, cf. Fig. 5a). During epoch I, frontal thinning rates were similarly intense at Kongsvegen and Kronebreen, but slightly less pronounced on Kongsbreen (Fig. 5b). Conversely, during epoch II, the largest elevation- change rates occur on Kongsvegen and Kongsbreen (about —2 to — 3ma^-1) while Kronebreen shows an elevation change of about —1 ma^-1 over the entire surface (Fig. 5c). Kronebreen’s frontal thinning rate in epoch II is only 33-66% of that observed in epoch I. The spatial and temporal ∂h/∂t variation between the glaciers and epochs suggests variation in ice dynamics.

Due to the lack of data in the upper area of Kronebreen (mainly in the 1990 DEM), the center-line differential GPS (dGPS) profile from 1996 is compared with both the 1966 and 2007 DEMs (Fig. 5b, c and g). Before differencing, the mean difference between the dGPS and overlapping 1990 DEM is estimated and removed from the 1996 heights to provide an estimated 1990 surface. This assumes that the relative geometry remained constant between 1990 and 1996 and all elevation differences are uniform within this period. This seems to be the case, at least for the longer 19662007 differences which remain relatively constant from about —0.75 to +0.25 m a^-1 above 700 m a.s.l. However, from 1966 to 1990, the elevation-change rates above 700 m a.s.l. were closer to zero while after 1990 they were closer to —1 m a^-1. Above 400 m a.s.l., the elevation change rate is less negative during epoch I than during epoch II, while below 400 m a.s.l. the elevation change rate is less negative during epoch II than in epoch I. Compared with Kongsvegen, there is less lateral variability of elevation changes at the highest elevation bins, suggesting that center-line values represent those for the entire elevation bin at Kronebreen but not Kongsvegen.

Surface mass-balance model

Table 3 shows calibrated DDF values and the rmse between measured and modelled SMB from various periods. On Kongsvegen, the 21 year stake record results in DDFs of 3.0 and 3.5 mm K^-1 a^-1 for snow and ice, respectively. On Kronebreen, DDFs for snow and ice are 3.1 and 4.7 mm K^-1 a^-1, respectively. Using these parameter sets, the model is applied to the time series of meteorological data. Figure 6 shows the mass balance measured at each stake vs the corresponding modelled value for all years (individual years are shown in the Appendix, Fig. 9). On Konsgvegen, the rmse of winter, summer and net stake mass balances are 0.24, 0.26 and 0.35 ma^-1, respectively (Table 3). On Kronebreen, the corresponding rmse values are 0.18, 0.29 and 0.33 m a^-1, respectively.

Fig. 6. The measured vs modelled winter (blue), summer (red) and net (black) stake mass balance on Kongsvegen from 1987 to 2008 and on Kronebreen from 2003 to 10. The rmse for each season is presented in Table 3.

Figure 7 shows the modelled cumulative net mass balance since 1969. In general, the SMB of Kongsvegen remained close to zero until becoming more negative in the late 1990s. Kronebreen was slightly positive during epoch I, but also turned more negative in epoch II. Superimposed on the figure are the geodetic balances estimated for the entire period and in epochs I (1966—1990/95) and II (1990/95—2007). The difference between the entire geodetic time series (1966— 2007) and the sum of epoch I and epoch II is the unremoved systematic bias shown in the residual triangulation of Table 2.

Fig. 7. The cumulative SMB model for Kongsvegen (red) and Kronebreen (black). The error zone is one standard deviation, σ, of model runs using DDFs of ±0.5mmK−1 d−1 around the central values provided in Table 3. The cumulative geodetic balance is also shown at squares and circles with connecting lines. The dashed red line shows the geodetic balance estimated using the full-range ∂h/∂t field, while the straight lines show those using the center line. Kronebreen is only shown with the center-line value as there was no visible difference between the center-line and full-range ∂h/∂t field for 1966–2007. All changes are relative to 1969, when the meteorological observations were initiated.

On Kongsvegen, the area-averaged elevation change, estimated by integrating ∂h/∂t over the glacier area, is m more negative than SMB (Fig. 7). The underestimation of the SMB model is an extrapolation failure, due to the predominance of z for distributing T and P, so horizontal variability is not accounted for. Spatial variability is seen in the elevation-difference maps (Fig. 5a and f) in which the center-line ∂h/∂t shows much more thickening than the average per elevation bin above the ELA. Integrating the center-line ∂h/∂t using the ‘hypsometric’ approach results in an area-averaged elevation change more similar to the modelled SMB (Fig. 7). Considering velocities are negligible on Kongsvegen (Reference Melvold and HagenMelvold and Hagen, 1998), the spatial variability of ∂h/∂t must be governed by SMB, most likely the distribution of snowfall, which is slightly positively biased towards the center line.

Figure 7 also presents the results of the SMB model on Kronebreen from 1969 to present. The cumulative mass balance remained slightly positive until the mid-1990s, when, similarly to Kongsvegen, the mass balance became more negative. The geodetic balance of Kronebreen shows a loss of ~18—20mw.e. within the time period. We suggest the strikingly large difference between the geodetic balance and SMB is caused by calving.

Flux divergence

Figure 7 makes it clear that on Kronebreen the residual between the geodetic balance and SMB is significantly negative. This component reflects the loss of ice due to calving, and is on the order of ten times the cumulative SMB. To analyze this in more detail, the average annual SMB and elevation change is plotted against elevation in Figure 8. On Kongsvegen, the difference between these curves is small and represents any systematic error related to ∂h/∂t and b. This error can be due, for example, to a failed prescription of p _eff in κ (Eqns (6) and (7)) or to unaccounted internal accumulation. On Kronebreen, the difference between ∂h/∂t and b is greater. We propose that is dominated by the emergence velocity below the ELA, as errors related to p _eff are minimal. Above the ELA, this may represent the submergence velocity, though greater uncertainty prevails due to the compressible firn area. The slope of with elevation increases from epoch I to epoch II (Fig. 8). This is interpreted as a larger influx of ice that compensates for ice loss through the SMB, resulting in a relatively flat ∂h/∂t gradient with elevation. Above the ELA, submergence rates may be on the order of —1 to —2 m a^-1, though unaccounted internal accumulation in the SMB or a lower k in the conversion of ∂h/∂t to water equivalent will reduce this value.

Fig. 8. The annual average surface mass-balance rates, b (red), elevation change rates, ∂h/∂t (black), and the difference between them (blue) on Kongsvegen (top) and Kronebreen (bottom). Elevation change rate pixels are shown in gray. On Kongsvegen, is essentially zero and thus the blue line represents an error term as ∂h/∂t and b have been shown to be equal (Reference Melvold and HagenMelvold and Hagen, 1998; Reference Hagen, Eiken, Kohler and MelvoldHagen and others, 2005). On Kronebreen, is positive below the ELA and negative above. The slope of with elevation increases from epoch I to epoch II.

The use of the smallest glacier area as a ‘reference’ surface for integration of the geometric changes and mass-balance components makes a drainage basin balance assessment possible (similar to Reference BindschadlerBindschadler, 1984). The 2007 front position functions as a flux gate in all calculations. The integration of the SMB over this area in both epochs is the glacier-wide balance (B in Eqn (13)). Integration of the elevation changes over this same area results in the volume change, ∂V/∂t. The difference results in an estimate of the ice flux, Q, through the 2007 flux gate. To account for the retreat, we also estimated the volume of ice loss over the area of retreat, ∂V _r/∂t, using the bathymetry in the retreat areas (Reference Kehrl, Hawley, Powell and Brigham-GretteKehrl and others, 2011). SMB of the retreat area, B _r, is estimated by applying the SMB model to the entire area of retreat, and linearly scaling the time series of B _r down to zero. This essentially assumes a linear retreat of the front.

Table 4 provides estimates of the components in Eqns (1214) for both glaciers and epochs assuming κ = 0.9. On Kongsvegen, B is the largest component of ∂V/∂t and the magnitude of B increased threefold from epoch I to epoch II. A slight residual between ∂V/∂t and B persists on Kongsvegen. Kronebreen exhibits the opposite behavior. In epoch I, B is 2.5% of ∂V/∂t, making the calving flux, Q, 97.5% of ice loss during this epoch. In epoch II, B doubled, comprising 25% of ∂V/∂t while Q represents 75% of ice loss. Despite the decrease in the proportion of calving to volume change, the effective calving through the 2007 flux gate, Q, increased by 1.5 times in epoch II compared with epoch I. The difference is slightly reduced when calculating the net flux, Q′, including the retreat flux, Q _r, which is dependent not only on Q but also on the underlying topography of the retreat area. Finally, Q′ and Q _r also include any ocean-induced melt that occurs on the vertical front of the glacier below sea level.

Assumptions and errors

The significance of errors on estimations is strongly dependent on the magnitude of the observed glacier changes. Many approaches for estimating errors in both in situ glaciologic mass balances and geodetic volume change estimates have been applied, commonly using the statistical theory of stochastic error propagation. For example, Reference Thibert, Blanc, Vincent and EckertThibert and others (2008) and Reference ZempZemp and others (2010) keep very detailed track of systematic uncertainties that are combined with the stochastic error to form the total error. The stochastic errors of the estimates are relatively easy to handle, while the systematic errors are difficult to detect and quantify. In this section, first we describe the stochastic and controlled systematic errors of each of the estimated mass continuity components. Then we discuss the systematic errors related to the assumptions we have applied, including sensitivity tests on these assumptions.

For the geodetic estimates, the stochastic point error is estimated as the standard deviation, σ, of elevation residuals between two DEMs over stable terrain (Table 2, maximum 10m). The error about the mean elevation change is then obtained by the standard error equation, assuming an autocorrelation length of 1 km (Reference Nuth, Kohler, Aas, Brandt and HagenNuth and others, 2007). Systematic errors in the DEMs are controlled by first coregistering the data, and the triangulated elevation residual between three DEMs over stable terrain is an estimate of the systematic bias remaining (the lower part of Table 2, specifically 1960-1990-2007 for Kronebreen and 19661995-2007 for Kongsvegen). Systematic and stochastic errors are combined through the root sum of squares (rss), resulting in a final error dominated by the systematic component. The geodetic error is then ~10-25% of the volume change (Table 4), though it can be >100% when the change is small (Kongsvegen, epoch I).

Errors of the in situ stake measurements at an annual resolution are relatively small compared with the error introduced by spatial extrapolation (i.e. not capturing the horizontal variability). Additionally, we model SMB empirically by relating the specific seasonal winter and summer measurements to temperature and precipitation. The accuracy of this approach is strongly dependent upon the calibrated model parameters. Therefore, error of the cumulative mass-balance series is estimated by taking one standard deviation of model outputs by varying DDFs for snow and ice of 0.5 mm °C^-1 a^-1 around the central value provided in Table 3. The error in B is ~30-200% of the Kongsvegen ∂V/∂t, since changes are not strongly different from zero. On Kronebreen, error in B is on the order of 1015% of ∂V/∂t (Table 4).

The error for the calving flux is then the combined rss of the ∂V/∂t and B errors, and is ~25% of the estimated calving flux on Kronebreen (Table 4). On Kongsvegen, we do not expect much ice loss through calving, so any residual between ∂V/∂t and B is associated with some unremoved systematic biases. Some ice could be lost through entrainment into Kronebreen at the tongue, though the sign of Q would be negative rather than positive. Interestingly, the magnitude of the bias is of the same order of magnitude as the estimated decadal average accumulation of superimposed ice (Reference Brandt, Kohler and LuthjeBrandt and others, 2008), an effect that would cause a positive value in the difference between ∂V/∂t and B.

The weakest part of our application of the mass continuity equation is the potential bias induced by our simplifying assumptions. The density ratio, k, required to convert elevation changes into water-equivalent volume changes is not completely certain and commonly assumes a constant value equal to the density of ice. This assumption is satisfied in the ablation area, but is questionable in the compressible firn pack. In the firn area, the density of ice can be assumed, based upon Sorge’s law, which states that the density at a given depth does not change with time, provided there is a constant accumulation rate (Reference BaderBader, 1954). However, a change in the firn-pack thickness and/or density will invalidate Sorge’s law. We also assume that englacial and basal mass-balance components are negligible, but these components are probably orders of magnitude smaller than the surface component (Reference Cuffey and PatersonCuffey and Paterson, 2010). Finally, mass-balance measurements, b, made in the accumulation area assume that all melt discharges from the glacier, although internal accumulation may be occurring. To assess the sensitivity of these assumptions on Q, we apply two scenarios, the first related to ∂h/∂t in the definition of k and the second related to b. The scenarios are based upon varying these terms above the average ELA, defined at 500 and 700ma.s.l. for Kongsvegen and Kronebreen, respectively.

• Scenario 1: Test of firn-thickness changes. Use k = 0.55 for dh/dt above the ELA.

• Scenario 2: Test of internal accumulation. Assume that 50% of the modelled mass losses, b, above the ELA are maintained in the system.

These scenarios provide upper and lower bounds to the derived flux through a fixed gate, Q, based upon simplified assumptions of systematic bias (Table 5). The effect of scenario 1 is opposite on the two glaciers because Kongsvegen is thickening above the ELA whereas Kronebreen is thinning above the ELA. Scenario 2 obviously increases the SMB. Cases of firn compaction would theoretically reduce k with limits approaching zero if the entire ∂h/∂t is the result of firn compaction. However, this is rather unlikely, mainly because the positive and negative observed ∂h/∂t above the ELA on Kongsvegen and Kronebreen, respectively, would imply that the firn pack of these adjacent glaciers is expanding and compressing, respectively, despite the similar climatic situations. In addition, no changes in the firn thickness of the upper basin of Kronebreen were observed between 1992 and 2005 (Uchida and others, 1993; Reference SjogrenSjogren and others, 2007; Reference Nuth, Moholdt, Kohler, Hagen and KääbNuth and others, 2010). In summary, the worst-case scenario of these simplifying assumptions results in a 20% difference in estimates of Q for Kronebreen. On Kongsvegen, the effects are slightly larger, due to the smaller magnitude of ∂V/∂t; however, the good coherence between ∂V/∂t and B through Eqn (11) (Fig. 7) suggests a proper assumption of k.

Interpretation

The comparison between SMB and geometric changes is based upon mass continuity. We also utilize concepts derived by Reference Elsberg, Harrison, Echelmeyer and KrimmelElsberg and others (2001) about the difference between the ‘conventional’ mass balance, which is the actual mass change, and the ‘reference’ mass balance, which is integration over a constant geometry. In their study, the reference surface is held constant at the largest area, and positive corrections are applied, based upon the geodetic volume changes, since the glacier is thinning and retreating. In our study, the ‘reference’ mass balance is obtained through integration over the newest (i.e. smallest) area, and transformation to conventional mass balances is in the negative direction to account for frontal retreat, B _r. The use of the smallest area additionally allows for full balance assessment (as in Reference BindschadlerBindschadler, 1984) of the glacier basin, by using the most recent front position as a flux gate.

Kongsvegen provides control on the SMB model because it is in its quiescent phase of a surge cycle (Reference Melvold and HagenMelvold and Hagen, 1998; Reference Hagen, Eiken, Kohler and MelvoldHagen and others, 2005). However, Figure 5f shows the center-line elevation change rates on Kongsvegen are not perfectly representative of the entire surface of Kongsvegen. Figure 7 shows the cumulative geodetic balances computed by the ‘grid’ approach and by the ‘hypsometric’ approach using the center line. The cumulative surface mass balance matches closest to the center-line geodetic estimates. The lack of significant ice flow suggests that the spatial ∂h/∂t variability, which is most dominant above the ELA, is governed by SMB variability. This may be due to radiation variability, which can affect melt (Reference Arnold, Rees, Hodson and KohlerArnold and others, 2006), previously hypothesized to cause geodetic balances to be more negative than either traditionally measured mass balances or center-line extrapolated estimates on a nearby small valley glacier, Midtre Lovenbreen (Reference Rees and ArnoldRees and Arnold, 2007; Reference Barrand, James and MurrayBarrand and others, 2010). If this were so, then we could expect that the reduced radiation along the edges of the glacier due to topographic shading (Reference Arnold, Rees, Hodson and KohlerArnold and others, 2006) would reduce melt and thus a center-line estimate would be more negative than a glacier-wide estimate. The opposite is seen here, and in previous studies.

We suggest instead that lateral accumulation variability is responsible for the underestimation of the modelled SMB and center-line geodetic balance (Fig. 7). Centerline ∂h/∂t above 500 m on Kongsvegen is more positive than the elevation bin average (Fig. 5), and slightly larger accumulation is measured along the main center line than in the two northern basins measured within the past 10 years (J. Kohler, unpublished data). Accumulation affects the timing of albedo decrease within the summer mass-balance season and the DDF transition from snow to ice within the model. This results in the model underestimating melt in areas where snow accumulation is overestimated. Therefore, the combined effect of an overestimated accumulation and an underestimated melt results in an underestimated cumulative SMB (and center-line geodetic balance) as compared to the full spatially integrated geodetic balance (Fig. 7).

The SMB of Kongsvegen has decreased drastically during epoch II (Table 4), especially since the late 1990s-early 2000s (Fig. 7). On Kronebreen, a similar decrease in the SMB is modelled, though this is not as great because of the higher elevation range of the Kronebreen basin. Similar rapid (geodetic) mass-balance decreases have been observed on two small land-terminating glaciers in western and central Spitsbergen (Reference KohlerKohler and others, 2007). Since 2007, it seems that the cumulative SMB has become less negative on both glaciers, similar to other measurements of surface mass balance and elevation changes in Svalbard (Reference Moholdt, Hagen, Eiken and SchulerMoholdt and others, 2010a,b).

On Kronebreen, ∂V/∂t estimated from the center-line ∂h/∂t and full spatially integrated ∂h/∂t are similar within the errors. As ∂V/∂t is much larger than B, the flux of ice through the 2007 flux gate, Q, is the main contributor to ∂V/∂t over the past four decades. Interestingly, despite the increased proportion of B to ∂V/∂t from epoch I to epoch II (2.5% to 25%), Q has nearly doubled (Table 4). Including the calving retreat flux, Q _r, in the estimate for the total calving flux, Q′, slightly reduces the difference. However, Q _r, and therefore Q′, is dependent upon both the influx of ice upstream, Q, and the underlying topography. For example, Reference Vieli, Funk and BlatterVieli and others (2001) model a rapid retreat over a basal depression and suggest that the rapid retreat is more dominantly an effect of bed topography rather than changes in climate. Therefore, it may be more feasible in this study to compare values for Q out of a fixed flux gate. It remains unclear how the retreat flux, Q _r, affects Q during the retreat of epoch I.

We can only speculate on reasons for the calving flux increase, as in situ measurements of velocity covering the entire period are not available. Available measurements do show interannual variability (Reference Kääb, Lefauconnier and MelvoldKääb and others, 2005), although no drastic increases have been detected over the past four decades (Reference VoigtVoigt, 1966; Reference Lefauconnier, Hagen and RudantLefauconnier and others, 1994b; Reference Kääb, Lefauconnier and MelvoldKääb and others, 2005). The flux of ice comprises both deformation and sliding components. Deformation velocity reacts slowly, due to its dependence upon ice thickness and slope. Conversely, glacier sliding reacts faster in relation to the subglacial drainage conditions (Reference Iken and BindschadlerIken and Bindschadler, 1986; Reference KambKamb, 1987). Recently, accelerated glacier motion has been observed to be related to water pulses that exceed the drainage system capacity (Reference Bartholomaus, Anderson and AndersonBartholomaus and others, 2008). Reference SchoofSchoof (2010) suggests that increased sliding velocities are due to more frequent water pulses, rather than an increase in the melt volume. Our calving estimates are averages over decades and may be indirectly related to the temporal average decrease in B. Since Kronebreen is losing volume over the entire glacier basin, especially in epoch II, ice deformation is supposedly decreasing. Therefore, it can be suggested that the increase in calving flux is related to the sliding velocity, and possibly an increase in the decadal average frequency and magnitude of water pulses.