Simulating the climatic mass balance of Svalbard glaciers from 2003 to 2013 with a high-resolution coupled atmosphere-glacier model

Introduction Conclusions References


Introduction
The Svalbard archipelago has a glacierized area of ca.34 000 km 2 (Nuth et al., 2013), representing ∼ 4 % of the world's land-ice mass outside the Greenland and Antarctic Figures ice sheets.If completely melted, the glaciers on Svalbard could potentially contribute to sea level rise of 17 ± 2 mm sea level equivalent (SLE; Martin-Espanol et al., 2015).The archipelago has already experienced significant warming during the 20th century (Førland et al., 2011) and, with the expected retreat of the sea ice margin, further warming as well as precipitation increase is expected (Day et al., 2012).Projections presented in the latest assessment report of the IPCC (AR5) shows that annual-mean temperatures in this region could rise between 7 and 11 • C by the end of the 21st century under the RPC8.5 scenario, accompanied by a projected precipitation increase between 20-50 % (IPCC 2013).Svalbard glaciers are therefore expected to undergo significant changes during this century (Day et al., 2012;Lang et al., 2015a).However, reliable estimates of future glacier changes require that we are able to reproduce recent observations.Current model estimates based on global climate datasets (Marzeion et al., 2012(Marzeion et al., , 2015) ) show significantly more negative mass balance in this region than satellite altimetry and satellite gravimetry over the last decade (Moholdt et al., 2010;Matsuo and Heki, 2013).
Regional model estimates of surface or climatic mass balance of Svalbard glaciers have so far mainly focused on individual or a few glaciers, and typically been based on empirical or statistical models (Lang et al., 2015b).A number of dynamical downscaling simulations focusing on Svalbard glaciers also exist (Day et al., 2012;Claremar et al., 2012;Lang et al., 2015a, b).However, only two of these studies compare their output with mass balance observations: Day et al. (2012;hereafter DA12) compares precipitation from HadRM3 RCM (25 km) with surface mass balance estimates from Pinglot et al. (1999), and Lang et al. (2015b;hereafter LA15b) compares output from the MAR model (10 km) to Pinglot et al. (1999Pinglot et al. ( , 2001) ) as well as a number of altimetry and gravimetry studies of Svalbard glacier mass balance (Wouters et al., 2008;Moholdt et al., 2010;Nuth et al., 2010;Mémin et al., 2011).Both studies show fair agreement with multi-year accumulation records from Pinglot et al. (1999Pinglot et al. ( , 2001)).LA15b also finds Svalbard mean elevation changes in good agreement with satellite estimates, even though the differences are substantial for some regions.However, DA12 and LA15b Figures do not validate mass balance estimates on time scales shorter than 4 years, nor do they validate on spatial scales that can capture variations on individual glaciers.DA12 also suggests that a grid spacing of 1-5 km may be needed to simulate surface mass balance in the complex terrain that is typical for Svalbard.
In this study, we aim to further reduce the spatial and temporal gap between observations and models, as well as to utilize a more extensive set of in situ measurements available in this region for model validation.We apply a coupled atmosphereglacier mass balance model to the entire Svalbard region with a horizontal grid spacing of 3 km, thereby capturing both regional averages for the period 2003-2013 as well as temporal and spatial variations of individual glaciers.The results are validated with (i) observations from weather stations, (ii) mass-balance stakes from four glaciers, (iii) snow accumulation across Austfonna, measured by ground-penetrating radar (GPR), and (iv) satellite altimetry.To examine the importance of model resolution in this region we also compare results from domains with 9 and 3 km grid spacing, and for a selected month, precipitation results from 9, 3 and 1 km grid spacing domains.
Through this extensive model evaluation we aim at enhanced reliability of the simulated climatic mass balance (CMB) and related variables required to quantify the state of Svalbard glaciers.

Methods
In the following sections, we describe the two components of the coupled modelling system: the Weather Research and Forecasting model (WRF; Sect.2.1) and the climatic mass balance model (Sect.2.2), including optimizations made in this study for high-Arctic conditions.In Sect.2.3, we describe the different validation data and sites before clarifying comparison methods in Sect.2.4.Introduction

Conclusions References
Tables Figures

Back Close
Full The WRF model is a state-of-the-art mesoscale atmospheric model (Skamarock and Klemp, 2008) widely used for research and forecasting applications.In Svalbard, the model has been used to study both atmospheric boundary layer processes (Kilpeläinen et al., 2011;Kilpeläinen et al., 2012), and atmosphere-land surface interactions over both tundra (Aas et al., 2015) and glaciers (Claremar et al., 2012), with horizontal grid spacing ranging from sub-kilometer scales to several tens of kilometers.In this study, we use the advanced research WRF version 3.6.1 configured with two nested domains of 9 and 3 km horizontal grid spacing.For a single month we also simulate the main regions of interest with additional nested domains at 1 km.The WRF model setup and forcing strategy follows that of Aas et al. (2015).The outer domain (9 km) covers a region of 1080 km × 1080 km using ERA-Interim (Dee et al., 2011) as lateral boundary conditions.It is one-way nested down to the 3 km domain covering all of Svalbard (Fig. 1).Within both domains the model is allowed to freely evolve (i.e.no nudging or re-initialization), and sea surface temperatures and sea ice fractions are prescribed based on the OSTIA dataset (Donlon et al., 2012).The physical parameterization options in WRF follow Aas et al. (2015) with the exceptions of the boundary layer and surface layer parameterizations, the vertical model resolution, and the use of explicit 6th order horizontal advection diffusion, which are selected following Collier et al. (2013).In addition, we use the newer NoahMP land surface scheme to simulate surface conditions and fluxes at non-glaciated grid cells, as it includes improved snow physics and multiple layers in the snowpack over the original Noah scheme (Niu et al., 2011).
To reduce computational time, we perform two five-year, transient simulations (2003-2008 and 2008-2013), both with a one-day spin-up for the atmosphere.We initialize glacier grid cells using results from an earlier simulation with a non-optimized version of WRF-CMB covering the same period.The first period has been initialized with the results for the year 2007, as no spin-up results are available for 2003, whereas the Introduction

Conclusions References
Tables Figures

Back Close
Full second period has been initialized with results for that year (2008).Comparing the initialization of the second period (starting in 2008) and the final results from the first period show there is good agreement for firn and snow outlines, making it far superior to the default uniform initialization.The one-month sensitivity simulation with 9, 3 and 1 km grid spacings has been performed as an isolated additional simulation with initialization directly from ERA-Interim to avoid problems with different spin up times for the different domains.

The climatic mass balance model
For glacier grid-cells, a modified version of the glacier climatic mass balance (CMB) model of Mölg et al. (2008Mölg et al. ( , 2009) ) is used to simulate glacier processes and surface fluxes.The CMB model computes the column specific mass balance from solid precipitation, surface and subsurface melt, refreezing, and liquid water storage in the snowpack, and surface vapor fluxes.The model solves the surface energy balance to determine the energy available for surface melt, and resolves the glacier subsurface down to a depth of 20 m.For glacier grid cells, the CMB model updates surface fields used in the WRF model (e.g., temperature; surface roughness, and albedo).Further information about the interactive coupling between the CMB model and WRF is given by Collier et al. (2013Collier et al. ( , 2015)).We make two adjustments to the CMB model to improve its suitability for Svalbard conditions.First, the albedo scheme (originally based on Oerlemans and Knap, 1998) has been modified to include a separate value for firn albedo (in addition to the standard categories of fresh snow, old snow, and ice).Firn is defined here as snow remaining from the last summer.Secondly, we introduce different aging models for snow (which determines the transition from fresh to old snow) for melting and sub-melting temperatures, by multiplying the snow aging parameter by a factor (warmfact, Table 1) when skin temperatures are at the melting point.The various albedo parameters (Table 1) have been selected based on observations at Austfonna, and thereafter adjusted to Introduction

Conclusions References
Tables Figures

Back Close
Full improve the simulated summer mass balance compared to in situ observations during one test year (2005)(2006); see also Sect. 6.1.1).
Finally, we note that some glacier grid-cells are affected by a numerical issue causing unrealistic sub-surface melting in the bottom snow layer, leading to unrealistically low corresponding CMB values.We were not able to pinpoint and remove the exact problem and therefore excluded sub-surface melting during the months October to May.This greatly reduced the effect of this problem with very little effect on all other grid cells.When averaged over the entire glacier area, this correction gave a CMB contribution of 1-17 mm w.e.yr −1 , with an average value of 8 mm w.e.yr −1 (less than 1 % of the average annual surface melt).This can therefore be considered as a minor correction to the total CMB, even if it is locally important.

Validation data and sites
Table 2 provides an overview of the datasets used for model evaluation.We selected four glaciers, described in more detail below, all of which have mass balance stake measurements at six or more locations, covering all or most years in our study period.The selected sites represent different conditions in terms of glacier geometry, geographical location, local meteorology, altitudinal range, and spatial extent (Fig. 1).Annual measurements of stake heights above the snow surface, snow depth, and density yield specific values of summer and winter mass balance (b s and b w ), which are combined to give the specific net mass balance, b n , or surface mass balance (SMB) for each balance year (i.e. between two consecutive end-of-summers).GPR measurements of snow accumulation along several transects across the Austfonna ice cap were made each spring in the period 2004-2013.Snow water equivalent (SWE) values are derived from the radar estimated snow depths multiplied with snow density determined at several snow pits, as described in more detail by Dunse et al. (2009).
Meteorological records at hourly resolution are available from two AWSs, one at Austfonna and one on Kongsvegen.Both datasets contain some data gaps.We ex-5781 Introduction

Conclusions References
Tables Figures

Back Close
Full tract hourly measurements of air temperature (∼ 2 m) and radiation (incoming and outgoing short-and long-wave), to calculate daily means, excluding days with incomplete records.Albedo is calculated as the ratio of the outgoing (SW out ) to incoming (SW in ) shortwave radiation, excluding observations outside the range [0.15, 0.95] or with SW in < 10 W m −2 .We estimate daily mean albedo using the five hourly observations closest to solar noon, to minimize the effect of low solar angle with associated large variations in measured albedo (Schuler et al., 2013).
A geodetic mass-balance estimate from repeat satellite altimetry for the period 2003-2008 (Moholdt et al., 2010) serves as an independent validation of the surface height changes due to climatic mass balance processes.However, it should be kept in mind that the geodetic mass-balance reflects both climatic and dynamic mass balance, i.e. it includes mass transfer from higher to lower elevations and losses due to calving at marine termini.We use the regional mean values between 2003 and 2008 according to Fig. 1.

Austfonna ice cap, northeast Svalbard
The Austfonna ice cap in the northeastern part of Svalbard (centered at 79.7 • N, 24.0 • E) is the largest ice cap of the archipelago.It covers an area of 7800 km 2 and has a simple dome-shaped topography, rising from sea level up to an elevation of ∼ 800 m a.s.l.(Moholdt and Kääb 2012).The recent CMB of Austfonna was nearly in balance (Moholdt et al., 2010), yet the ice cap was losing mass due to calving and retreat of the marine margin (Dowdeswell et al., 2008).Snow accumulation is spatially and temporally heterogeneous; accumulation is asymmetrical across the ice cap, with amounts in the southeast being double the amounts in the northwest, and there is large interannual variability along all profiles (Pinglot et al., 1999;Taurisano et al., 2007;Dunse et al., 2009).
Since spring 2004, field measurements have been performed annually by the University of Oslo and the Norwegian Polar Institute.Available data include records from about 20 mass balance stakes, annually repeated GPR and kinematic GNSS (Global Naviga-5782 Introduction

Conclusions References
Tables Figures

Back Close
Full tion Satellite System) profiling across the ice cap, and snow pit investigations of snow depth and density (Taurisano et al., 2007;Dunse et al., 2009).In the present study we compare GPR-derived winter accumulation with the corresponding WRF-CMB results, averaging all available in situ measurements within a particular WRF-CMB grid cell (Sect.3.3).
Etonbreen is located at the western part of Austfonna, with six stakes, and an AWS operated since 2004.The AWS is located at 22 • 25 12 E, 79 • 43 48 N and 370 m a.s.l., just below the mean equilibrium line at ∼ 400 m (Schuler et al., 2013).
Kongsvegen is a ∼ 100 km 2 , ∼ 27 km long valley glacier extending from an ice divide at ∼ 800 m a.s.l. down to sea level.Outflow at its marine terminus is restricted by its fast-flowing neighbor Kronebreen, with which Kongsvegen shares a small fraction of the calving front.The Norwegian Polar Institute has measured winter and summer mass balance at nine stakes since 1986 (Hagen et al., 2003;Nuth et al., 2012;Karner et al., 2013).Kongsvegen is a surge-type glacier, currently in its quiescent phase, since the last surge around 1948.Observed elevation changes are dominated by the SMB (Melvold and Hagen, 1998).
North of Kongsvegen is Holtedahlfonna, the upper catchment of the Holtedahlfonna-Kronebreen glacier system, whose total area is ∼ 390 km 2 , and which extends up to an elevation of ∼ 1400 m a.s.l.SMB has been studied on Holtedahlfonna since spring 2003, using ten mass-balance stakes.
Stakes at both glaciers are measured in nearly all years in spring and at the end-ofsummer, typically in early September.During the last two decades, the SMB of both glaciers has changed from close to zero to increasingly negative values (Kohler et al., 2007;Nuth et al., 2012).Since spring 2000, an AWS has been operated on Kongsve-Introduction

Conclusions References
Tables Figures

Hansbreen, southern Spitsbergen
Hansbreen (77.1 • N, 15.6 • E) is located in the southern part of Spitsbergen, and covers an area of ∼ 56 km 2 .It is a 15 km long valley glacier, extending from its ∼ 1.5 km wide active calving front in Isbjornhamna, Hornsund, up to an ice divide at ∼ 490 m a.s.l.
The Hansbreen system consists of the main trunk glacier and 4 tributary glaciers on the west side (Grabiec et al., 2011).SMB on Hansbreen is studied since 1989 when 11 mass balance stakes were deployed along centerline of the glacier.Since 2005 these stakes are measured on weekly basis in ablation zone and on monthly basis in accumulation area (Grabiec et al., 2012).Snow accumulation on Hansbreen is highly variable.Earlier studies show that there is a strong asymmetry between eastern and western side of the glacier.Minimal snow accumulation is observed in southern part of the tongue and along eastern side of the glacier.On the east side, the glacier is bordered by the massif of Sofiekammen, which forms an orographic barrier for advecting air masses.Therefore a foehn effect is widely observed during strong easterly winds, causing deflation of snow and redeposition towards the western side of the glacier (Grabiec et al., 2006).

Comparison methods
We compare our model results with AWS (Sect.3.1) and stake (Sect.3.2) data using the WRF-CMB grid point nearest to each data point.All stakes on a glacier are compared with the corresponding grid cells to yield information both about mean values as well gradients along the glacier.The ELA (also Sect.also includes internal refreezing below the last summer surface (LSS), which is not captured by the stake SMB.

Model validation
The main goal with this study is to evaluate the ability of the WRF-CMB model to reproduce measured CMB and height changes at Svalbard over the study period, both on the regional and local scales.In the following sections, we assess the model performance by comparing the simulated SEB and temperature with the AWS data (Sect.3.1); the CMB at the four glaciers with stake measurements (Sect.3.2); the winter accumulation at Austfonna with GPR measurements (Sect.3.3); and, finally, the simulated regional height changes over the first five years of the study period with the altimetry data (Sect.3.4).All model results in this section are from the 3 km domain.

Weather stations
The AWSs provide key information about the conditions at the glacier surface throughout the year, and therefore permit a detailed evaluation of important aspects of the model.Table 3 compares  A more detailed comparison of simulated and observed albedo from the summers of 2008 and 2013 is shown in Fig. 2, along with simulated solid and liquid precipitation.Overall, the model simulates well both the magnitude and temporal changes in albedo.The close connection between simulated solid precipitation events and observed albedo increase indicates that the model adequately captures both the timing and phase of the summer precipitation.The simulated albedo response to solid precipitation is also similar to measurements, except for at Etonbreen in 2013.Here the model does not simulate a large enough increase in albedo after snow events, which might indicate that the model is not sensitive enough to small amounts of snow on ice (i.e.too large depth scale in Table 1), or simulates too little snow on these occasions.
Additionally, the model underestimates snow albedo during much of the summer in 2008.While the observations show a very slow decrease in albedo over this summer, the modeled albedo quickly drops to values below 0.7.This could be related to the almost complete lack of rain events during the summer of 2008.In 2013, on the other hand, the model simulates numerous and large rain events during the summer, which seems to coincide with observed drops in albedo.This suggests that snow wetness needs to be accounted for in the albedo parameterization to realistically simulate snow albedo on Svalbard, and that the snow aging parameters used here (Table 1) are more appropriate for wet rather than dry summer conditions.
Altogether, the model seems to reasonably reproduce the local conditions at the AWS stations.Individual radiation biases are found, which are likely to impact the quality of the CMB estimations negatively.However, this is a known challenge in the Arctic, where atmospheric models often show significant biases in key cloud properties (Morrison et al., 2009).as having anomalously low and high mass balance, respectively (Fig. 3, stippled lines).The agreement, however, is not as good as for the mean values.

Stakes
To look further at temporal variations, we compare modeled ELA to that derived from stake measurements (Fig. 4).In general there is good agreement both in terms of ELA height as well as inter-annual variability, with the exception again being Hansbreen.Here the model overestimates ELA, in accordance with the underestimation of b w (Fig. 3).At Hansbreen glacier the observed mass balance -elevation relationship reveals considerable non-linearity, and the simple ELA estimation (Sect.2.4) used here might be inappropriate.For the other three glaciers the model simulates ELA well, including the large difference between Kongsvegen and Holtedahlfonna, which are located close to each other (Fig. 1).

GPR-derived snow profiles
To evaluate the winter-mass balance simulated by WRF-CMB we compare simulated precipitation across Austfonna between early September and early May with observations of snow accumulation by GPR.The results from May 2004 and May 2006 (Fig. 5) demonstrate the large inter-annual variability in the total amount of snow, with 2004 and 2006 representing a low (Fig. 5a) and high (Fig. 5b) accumulation year, respectively.WRF-CMB captures the spatial pattern and total amount of snow very well, with average biases of less than 1 and ∼ 3 % and mean absolute errors (MAEs) of 0.05 and 0.06 m w.e.yr −1 in 2004 and 2006 respectively.Local differences, for example an overestimation of snow within the ablation area towards the lower end of the western profile in both years (Fig. 5b and d accumulation in the summit area, which again could be related to wind redistribution of snow, or by a negative elevation bias of the WRF-CMB DEM which is on average ∼ 20-50 m lower than the average position of the GPR measurements.

Satellite altimetry
We perform a regional evaluation of the simulations by comparing surface height changes with those measured by satellite altimetry in the period 2003-2008 (Moholdt et al., 2010).Note here that measured height changes include the effect of glacier dynamics, which are not simulated by the model.Still, the comparison shows good agreement for each of the six regions, as well as for Svalbard as a whole (Table 4).Both model and satellite data show Austfonna and northeast Spitsbergen to be the only regions with positive surface height change during these years, and northwest Spitsbergen as the region with the largest surface lowering.The other three regions all show moderate lowering in both estimates.Together, these results show that the model captures both the mean CMB value, as well as its spatial and temporal variability across Svalbard well.The largest model errors are found at Hansbreen, but this glacier also has the largest cross-glacier variability in accumulation, as well as a relatively steep accumulation gradient, such that the stake measurements may not correlate as well to the model as at other sites with a greater degree of spatial homogeneity.This raises the question about the sensitivity of model resolution in relation to capturing the main topographic features.

Sensitivity to model resolution
Svalbard topography is relatively rugged, with fjords, tundra, glaciers and mountains all found in close proximity to each other.Dynamical downscaling requires a discrete representation of the topography and has therefore limited spatial resolution, which potentially can be insufficient to resolve a number of small-scale processes.In the fol-Introduction

Conclusions References
Tables Figures

Back Close
Full lowing section, we investigate the sensitivity to model resolution by comparing CMB from the 9 and 3 km domains for the entire model period.We then evaluate precipitation amounts and distribution in the simulations with 9, 3 and 1 km grid spacings (hereafter R9, R3 and R1, respectively) in the selected month of October 2007, which was among the wettest months in the 10-year period, with a number of smaller and larger precipitation events.
Figure 6 shows the mean annual R9 and R3 CMBs over the entire period.When averaged across all of Svalbard, the difference in CMB is small (less than 10 mm w.e.yr −1 ).On a regional scale, however, differences are more substantial, with the southeastern islands (BE region in Fig. 1) having more than 100 mm w.e.yr −1 more negative mass balance in R9 than in R3.The small difference found when averaging across the entire archipelago therefore seems to be at least partly due to compensating regional differences.On the local scale, differences are even more pronounced.The lack of small glaciers in R9 gives large ice-free areas compared with R3 (e.g. in much of central Spitsbergen).The CMB gradient is also often larger in R3 than R9, as higher maximum and lower minimum values are resolved.For long-term model simulations, where geometry adjustments of glaciers are to be considered, resolving the adequate range of CMB would be important.
The level of detail resolved also plays an important role in model evaluation, especially when comparing to in situ point measurements.Figure 7 shows the terrain height in the three regions with stake measurements, at the resolutions R9, R3 and R1.Here we see new topographic features and a more detailed coastline emerging with each increase in resolution.In R9 the topographic gradients along the glaciers with stakes are mostly too low, and many stakes maps on to the same grid cell.Conversely, in R3 most stakes maps on to individual grid cells and the altitudinal range is in better agreement with the stakes (see also Fig. 3).Still, there are distinct topographic features around these four glaciers that only emerge in R1.For example, the model represents Hansbreen, Kongsvegen and the upper part of Holtedahlfonna as valley glaciers only at this resolution.The simulated precipitation from October 2007 Introduction

Conclusions References
Tables Figures

Back Close
Full with the three resolutions (Fig. 8) also shows distinct differences.At Etonbreen, the shape of the precipitation curve along the glacier is very different in R9 compared with R3 and R1.The other three glaciers show more linear increases in precipitation with altitude during this month at all three resolutions.Still the precipitation amount varies and mainly increases with higher resolution, with local differences of up to about 50 mm (25 %), at Hansbreen (∼ 300 m a.s.l.).We therefore conclude that increasing resolution from 9 to 3 km grid spacings is important for simulating precipitation and CMB on local and regional scales, as well as for reliable comparison with in situ point measurements.Increasing the resolution further to 1 km grid spacing seems to have a smaller effect on the four glaciers investigated, although the precipitation amount is further increased.Model resolution might therefore explain some of the negative model bias in b w (cf.Fig. 3).However, resolution alone does not explain the large b w bias at Hansbreen; other processes are important for the accumulation pattern at this glacier, for instance wind drift and snow redistribution.

Climatic mass balance
We now turn to simulated annual and seasonal CMB results for Svalbard and the different sub-regions (Fig. 1).As has also been reported in other studies (e.g.Hagen et al., 2003;Lang et al., 2015b) we find large inter-annual variation in CMB (Fig. 9), which precludes trend analysis on the timescales considered in this study.Likewise, our mean mass balance value of −167 mm w.e.yr −1 is dependent on the model pe-  [2006][2007][2008][2009][2010].Multiplied with a glacier area of 34 000 km 2 (Nuth et al., 2013) this corresponds to a mean annual mass loss of 5.7, 7.7 Gt, or a mean annual mass gain of 1.6 Gt, respectively.
Comparing the two seasons reveals that summer mass balance varies about twice as much as winter mass balance (Fig. 9), indicating that ablation processes vary more from year to year than winter accumulation.This is confirmed from the annual mass Introduction

Conclusions References
Tables Figures

Back Close
Full fluxes shown in Fig. 10b, with solid precipitation showing much less variability than surface melting.Melting is in turn largely a result of the radiation imbalance during the summer months (JJA; Fig. 10a), which on average accounts for about 80 % of the net energy to the surface during these months.However, years with anomalously large melting (especially 2004 and 2013) are characterized by much larger than normal sensible and latent heat fluxes at the surface.As these fluxes require warmer and moister air in the atmosphere relative to the glacier surface, these melting anomalies likely result from advection of warm air from outside the region, as was also found by Lang et al. (2015b) for the year 2013.
For the individual regions we find large differences in simulated CMB that persists throughout the whole period (Fig. 9; right axis).Most noticeably the northeastern regions (AF, VF and NE) show less negative summer balance than the Svalbard average and mostly lower than average winter accumulation.The regions with high winter accumulation in the south and east correspond to the regions with the lowest ELA reported by Hagen et al. (2003; see also mean annual precipitation and ELA estimates in the Supplement).

Comparison with earlier studies
We will focus our comparison here with the two regional model studies that report comparison with surface or climatic mass balance measurements, namely DA12 and LA15b (Sect.1).Both of these studies compare their results with ice core measurements in accumulation areas (Pinglot et al., 1999(Pinglot et al., , 2001)).DA12 finds biases ranging from −0.24 to 0.13 m w.e.y −1 with a mean absolute error close to 0.1 m w.e.y −1 , whereas LA15b reports biases between −0.31 and 0.14 m w.e.y −1 , also with a MAE of 0.1 m w.e.y −1 .Our simulation period does not cover these ice core measurements so that a direct comparison is not possible.Still, the MAEs in winter accumulation at Austfonna (Fig. 5) are

TCD Introduction Conclusions References
Tables Figures

Back Close
Full at Hansbreen where it ranges from close to zero to about 1 m w.e.y −1 .Our results are therefore an improvement compared to DA12 and LA15b for CMB when Hansbreen is ignored (where processes not represented in our model are believed to be important).
The near surface temperature biases are in both DA12 and LA15b larger than ours.DA12 reports annual biases "less than 2 • C" for two stations, and "less than 6 • C" for the third (Kongsvegen) and LA15 showed annual biases between −1.3 and −4.0 • C (all negative).We only use data from the two AWSs on the glaciers here, with biases of 0.2 and −1.9 • C at Kongsvegen and Etonbreen respectively.These biases span the two-year mean biases found by Claremar et al. (2012), and are also similar to those found by Aas et al. (2015) for the year 2008-2009 (including also tundra sites), both using the WRF model.We note here that our results do not cover the same time periods as DA12 and LA15b, so that the quality of the boundary conditions might differ.Smaller elevation differences between the station and the corresponding model grid cell in our 3 km domain compared to DA12 and LA15b probably also contribute to the improvement seen here.However, LA15b reports relatively small biases for ERA-Interim (between −1.95 and 2.24 • C) despite similar or larger elevation differences and covering the same period.It therefore seems clear that the WRF-CMB model with the setup used here offers a real improvement over DA12 and LA15b.

Model uncertainties
In the following section we discuss the main uncertainties in the model results, starting with the atmospheric forcing, before discussing the representation of albedo, turbulent fluxes and sub-surface processes in the CMB model.Although we recognize that the observations have uncertainties and limitations, it is beyond the scope of this study to go into those here.

Conclusions References
Tables Figures

Back Close
Full

Atmospheric forcing
The quality of any CMB simulation depends largely on the atmospheric forcing used.In this study we have used a coupled atmosphere-glacier model in which the atmospheric component (WRF) has been well tested for this region.It has been shown to produce temperatures that are in good agreement with observations, and relatively small biases in energy fluxes on annual time scales, although they can differ significantly on seasonal and shorter time scales (Aas et al., 2015).Deviations from observations can come both from insufficient representation of processes within the WRF model and from errors in the boundary conditions (here ERA-Interim and OSTIA).However, comparison with weather stations (Sect.3.1) and the results from Aas et al. (2015) show that there is good agreement with observations for the atmospheric part, compared to similar studies (Claremar et al., 2012;DA12;LA15b).

Albedo
The albedo parameterization in the CMB model uses a set of parameters that cannot be universally set, but instead varies considerably for different locations.Initially, these parameters were set based on data from the AWS at Austfonna, where we had the longest series of radiation data.This, however, gave too little summer ablation in general, and these parameters were therefore instead tuned to give summer mass balance values in line with observations based on a set of sensitivity simulations of the year 2005-2006 with the 9 km grid spacing domain.The resulting values (Table 1) are similar to those used by van Pelt et al. (2012).With the complex model system used here we cannot perform a set of simulations at the full resolution and time period, but only acknowledge that these values are uncertain, and should ideally vary across the region.Greuell et al. (2007) found MODIS ice albedo values at different glaciers across Svalbard between 0.44 (Etonbreen) and 0.13 (Lisbethbreen), which also explains why our initial values from Etonbreen gave too little summer ablation in general.A large step forward for CMB modelling of this region would therefore be to include spatially varying albedo parameters based on satellite measurements.Including snow wetness could also offer improvements for simulation of changes in snow albedo with time (see Sect. 3.1).

Atmospheric stability
We have seen that the sensible and latent heat fluxes are an important part of the SEB during the ablation season (Sect.5).However, the simulation of these fluxes in stable conditions is a major challenge subject to ongoing research (e.g.Holtslag et al., 2013).
In our setup of the CMB model we use the bulk Richardson number to directly calculate the stability correction of the turbulent fluxes.Conway and Cullen (2013) found this to give too low fluxes during stable conditions and low wind speeds at a New Zealand Southern Alps glacier, with simulated turbulent fluxes close to zero when observations showed values between 50-100 W m −2 .Based on their results, and to avoid runaway cooling of the glacier surface during stable conditions in the winter, we had to limit the stability correction of turbulent fluxes to 30 % in stable conditions (Collier et al., 2015).Due to these issues, and since they are not measured at the study glaciers, these fluxes also contribute to the overall model uncertainty.

Sub-surface processes
As can be seen from the inclusion of the sub-surface processes likely offers an important step forward compared to only simulating the surface mass balance, the accuracy of these sub-surface processes in the model is uncertain.

Conclusions
In this study we have simulated the climatic glacier mass balance (CMB) of the Svalbard archipelago with a coupled atmosphere-glacier model, for the period 2003 to 2013 with 9 and 3 km grid spacings, as well as with 1 km grid spacing for a shorter period.
The results have been compared with extensive observational data from weather stations, mass balance stakes, ground penetrating radar and satellite altimetry.Our main findings are: -The WRF-CMB model with 3 km grid spacing and the configuration used here reproduces observed CMB on Svalbard very well, from local to regional scales.
-We confirm the need for high resolution (1-5 km grid spacing) to realistically simulate CMB at glacier scale, as suggested by Day et al. (2012).The results from 3 and 1 km grid spacings show distinct features at local scales that are not present at 9 km, and mean CMB at regional scale differed up to 0.1 m w.e.yr −1 between the 3 and 9 km simulations.
-Large variations in CMB on small spatial scales reduce the representativeness of individual point measurements when compared with grid cells larger than 1-5 km.
-We find large year-to-year variability in average Svalbard CMB during our simulation period, which can be mainly attributed to variations in melting.The largest component in the summer surface energy balance driving this melting is the radiation imbalance, even though the main feature separating the years with larger than normal melting is large sensible and latent heat fluxes.Introduction

Conclusions References
Tables Figures

Back Close
Full For the first time, this study presents results from dynamical downscaling of Svalbard CMB at the resolution suggested by Day et al. (2012).In addition, we have utilized a large number of observations on different spatial scales to get a robust evaluation of model performance, thereby representing a considerable step forward in the pursuit of reliable simulations of the CMB on Svalbard.Further improvements to several aspects of the WRF-CMB would be desirable for this region, including using spatially variable albedo parameters and improved representation of sub-surface processes.Still the WRF-CMB model -which has not previously been validated for Arctic conditions -has here been shown to be an appropriate tool for studying Svalbard CMB.
The Supplement related to this article is available online at doi:10.5194/tcd-9-5775-2015-supplement.Introduction

Conclusions References
Tables Figures

Back Close
Full  Full  Full Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | gen, at 78.76 • N, 13.16 • E at 537 m a.s.l., close to the equilibrium line altitude (ELA; 3.2) is calculated by linear interpolation of the two stakes or grid cells with the least positive and negative mass balance.When all stakes or grid cells have the same sign, we use the maximum and minimum value to extrapolate to the ELA.Note that the CMB simulated by WRF-CMB Introduction Discussion Paper | Discussion Paper | Discussion Paper | simulated and observed air temperature, radiation fluxes, and albedo at Kongsvegen and Etonbreen.At Kongsvegen the simulated temperature compares very well with observations, with a bias of less than 0.2 • C and correlation of daily mean values of 0.98.At Etonbreen the simulated temperatures are somewhat too low (−1.9 • C), but with a similar correlation (0.97).The radiation components also show better agreement at Kongsvegen than at Etonbreen.At Kongsvegen, biases ranging from 0.4 W m −2 (outgoing longwave, LW out ) to −6.6 W m −2 (incoming longwave, LW in ), whereas the radiation biases at Etonbreen vary from −4.0 (LW out ) to −12 W m −2 (SW in and LW in ).There is also a noticeable albedo bias of −0.08 at Etonbreen, compared to only −0.03 at Kongsvegen. Discussion Paper | Discussion Paper | Discussion Paper |

Figure 3
Figure 3 shows mean annual (b n ), winter (b w ), and summer (b s ) CMB at each stake over the study period.Overall, the model reproduces the mean observed mass balances well, including differences in mean values between the four glaciers and the 5786 ) may partly be explained by wind erosion, which is not represented in the model.The model also tends to slightly underestimate snow 5787 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | riod considered; five-year mean values vary from −227 mm w.e.yr −1 (2009-2013) and 48 mm w.e.yr −1 ( Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | 0.05 and 0.06 m w.e.y −1 for 2004 and 2006 respectively.At the stake locations the simulated CMB (Fig. 3) are mostly within about 0.1 m w.e.y −1 of the observations, except Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Fig. 10b the sub-surface components of the CMB (refreezing, superimposed ice, subsurface melting, and change in liquid water content) contribute considerably to the total simulated CMB.An important simplification in the CMB model is the use of a bulk snow density, as it is not intended for detailed snowpack studies.For the highest areas on Austfonna and Northeast Spitsbergen, where the simulated snow and firn depth reaches values of 14 m during the simulation, this simplification may be inappropriate.In addition, the CMB model currently calculates superimposed ice as a specified fraction of the internal refreezing of liquid water.Modeling this process is challenging, even with a vertically resolved treatment of snow density.Thus, while TCD Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Karner, F., Obleitner, F., Krismer, T., Kohler, J., and Greuell, W.: A decade of energy and mass balance investigations on the glacier Kongsvegen, Svalbard, J. Geophys.Res.-Atmos., 118, 3986-4000, 2013.Kilpeläinen, T., Vihma, T., and Olafsson, H.: Modelling of spatial variability and topographic effects over Arctic fjords in Svalbard, Tellus A, 63, 223-237, 2011Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper |

Figure 1 .
Figure 1.Main figure: land areas in the 3 km model domain.Colors indicate glacier grid cells in different sub-regions and gray indicate non-glacier land grid cells.Stake and AWS locations at the four main validation glaciers are shown as red * and black *, respectively.NW: northwestern Spitsbergen, NE: northeastern Spitsbergen, SS: southern Spitsbergen, BE: Barentsøya and Edgeøya, VF: Vestfonna and AF: Austfonna.Overlay figure: glacier hypsometry from 3 km model domain compared with 90 m DEM(Nuth et al., 2013).

Figure 2 .Figure 3 .
Figure 2. Simulated (blue) and observed (red) albedo at Etonbreen and Kongsvegen for the summers 2008 and 2013, with simulated solid (blue) and liquid (yellow) precipitation indicated as bars.

Table 2 .
Description of observations.

Table 3 .
Bias, mean absolute error (MAE) and correlation between simulated and observed daily mean temperature, radiation and albedo at Etonbreen and Kongsvegen.