Estimation of the Antarctic surface mass balance using the regional climate model MAR (1979–2015) and identification of dominant processes

Abstract

The Antarctic ice sheet mass balance is a major component of the sea level budget and results from the difference of two fluxes of a similar magnitude: ice flow discharging in the ocean and net snow accumulation on the ice sheet surface, i.e. the surface mass balance (SMB). Separately modelling ice dynamics and SMB is the only way to project future trends. In addition, mass balance studies frequently use regional climate models (RCMs) outputs as an alternative to observed fields because SMB observations are particularly scarce on the ice sheet. Here we evaluate new simulations of the polar RCM MAR forced by three reanalyses, ERA-Interim, JRA-55, and MERRA-2, for the period 1979–2015, and we compare MAR results to the last outputs of the RCM RACMO2 forced by ERA-Interim. We show that MAR and RACMO2 perform similarly well in simulating coast-to-plateau SMB gradients, and we find no significant differences in their simulated SMB when integrated over the ice sheet or its major basins. More importantly, we outline and quantify missing or underestimated processes in both RCMs. Along stake transects, we show that both models accumulate too much snow on crests, and not enough snow in valleys, as a result of drifting snow transport fluxes not included in MAR and probably underestimated in RACMO2 by a factor of 3. Our results tend to confirm that drifting snow transport and sublimation fluxes are much larger than previous model-based estimates and need to be better resolved and constrained in climate models. Sublimation of precipitating particles in low-level atmospheric layers is responsible for the significantly lower snowfall rates in MAR than in RACMO2 in katabatic channels at the ice sheet margins. Atmospheric sublimation in MAR represents 363 $\mathrm{Gt}$ ${\mathrm{yr}}^{-\mathrm{1}}$ over the grounded ice sheet for the year 2015, which is 16 % of the simulated snowfall loaded at the ground. This estimate is consistent with a recent study based on precipitation radar observations and is more than twice as much as simulated in RACMO2 because of different time residence of precipitating particles in the atmosphere. The remaining spatial differences in snowfall between MAR and RACMO2 are attributed to differences in advection of precipitation with snowfall particles being likely advected too far inland in MAR.

Introduction

Mass loss from the Antarctic ice sheet (AIS) and therewith its contribution to the sea level budget results from the difference of two fluxes of a similar magnitude: ice flow discharging in the ocean (D) and net snow accumulation on the ice sheet surface, i.e. the surface mass balance (SMB). The total ice sheet mass balance (SMB minus D) can be assessed using satellite altimetry, gravimetry, or the input–output method , which all request SMB estimates. The input–output method, which consists in separately modelling ice dynamics and SMB, is also the only way to project future trends.

SMB as used in this study is the sum of mass gains (mainly snowfall accumulation and some riming), mass losses (mainly surface and drifting snow sublimation, some liquid water run-off), and drifting snow transport (defined as the horizontal advection of the drifting snow), which can lead to either mass gain or mass loss. Snowfall rates are 1 order of magnitude larger than all of the other SMB fluxes at the continental scale , with the largest amounts found along the ice sheet margins due to cyclonic activity in the Southern Ocean and to the orographic lifting of relatively warm and moist air masses . Accumulation patterns are highly variable at the kilometre scale and from year to year e.g.. Consequently, proper observations of SMB require a high spatial coverage (e.g. stake lines, accumulation radars plus ice cores for layer dating and snow density) and a temporal sampling spanning several years . Even if efforts have been made to fulfil those requirements, ground-based observations are scarce and carry with them high logistical costs in this cold, windy, and remote environment. Interpolation techniques used to interpolate the scarce SMB observations encounter major caveats .

This is why many AIS mass balance studies use output of regional climate models (RCMs) to estimate ice sheet SMB for the recent decades e.g.. In order to obtain a good agreement with observations, atmospheric models require accurate large-scale circulation patterns together with a proper representation of snow surface processes, clouds, and turbulent fluxes and a relatively high horizontal resolution to properly resolve the complex ice sheet topography at the margins.

Here, we present new simulations of the RCM MAR, applied for the first time over the whole AIS, but already widely used for polar studies, e.g. in Greenland , Svalbard , Adélie Land Antarctic coastal area;, and Dome C Antarctic plateau;. We compare MAR-simulated SMB with the state-of-the-art RCM RACMO2 . We use available SMB observational datasets to show that MAR and RACMO2 perform similarly well in simulating the SMB spatial gradients. In addition, we identify significant processes that still need to be included or improved in both RCMs.

In Sect. , we describe MAR and its specific set-up for Antarctica, together with RACMO2, the forcing fields, observational datasets, and methods designed for model evaluation. In Sect. , we show that both RCMs share common biases against observed SMB, resulting from drifting snow transport fluxes. Secondly, we analyse SMB differences among models and show that many of the discrepancies can be attributed to low-level sublimation of precipitation in katabatic channels and to the difference in precipitation advection inland. Finally, in Sect. , we summarise our main findings and discuss further efforts to be achieved for a better assessment of the AIS SMB.

Data and methods

Regional modelling

Regional atmospheric models

For the first time, the polar-oriented regional atmospheric model MAR is applied for decades-long simulations over the whole AIS. MAR atmospheric dynamics are based on the hydrostatic approximation of the primitive equations, fully described in . Prognostic equations are used to depict five water species: specific humidity, cloud droplets and ice crystals, raindrops, and snow particles . Sublimation of airborne snow particles is a direct contribution to the heat and moisture budget of the atmospheric layer in which these particles are simulated. The radiative transfer through the atmosphere is parametrised as in , with snow particles affecting the atmospheric optical depth . The atmospheric component is coupled to the surface scheme SISVAT soil ice snow vegetation atmosphere transfer; dealing with the energy and mass exchanges among surface, snow, and atmosphere. The snow-ice part of SISVAT is based on the snow model CROCUS . It is a one-dimensional multilayered energy balance model which simulates meltwater refreezing, snow metamorphism, and snow surface albedo depending on snow properties. We used MAR version 3.6.4, simply called MAR hereafter. In this version the physical settings are the same as in MAR version 3.5.2 used for Greenland , except for the adaptations detailed below.

Grid. Projection is the standard Antarctic polar stereographic method (EPSG:3031). The horizontal resolution is 35 $\mathrm{km}$, an intermediate resolution that results from a computation time compromise in order to run the model with multiple reanalyses and global climate model forcings over the 20th and the 21st centuries. The vertical discretisation is composed of 23 hybrid levels from $\sim \mathrm{2}$ $\mathrm{m}$ to $\sim \mathrm{17}\mathrm{000}$ $\mathrm{m}$ above the ground.

Boundaries. The topography is derived from the Bedmap2 surface elevation dataset . Because the Antarctic domain is about 4 times larger than the Greenland domain, the circulation has to be more strongly constrained. This is why we use a boundary relaxation of temperature and wind in the upper atmosphere starting from 400 $\mathrm{hPa}$ ($\sim \mathrm{6000}$ $\mathrm{m}$ above the ground) to 50 $\mathrm{hPa}$ (upper level), as in , whereas relaxation starts from 200 $\mathrm{hPa}$ in .

Parameterisations.

• a.

  The surface snow density ${\mathit{\rho }}_{\mathrm{s}}$ ($\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{3}}$) is computed as a function of 10 $\mathrm{m}$ wind speed ws${}_{\mathrm{10}}$ (m ${\mathrm{s}}^{-\mathrm{1}}$) and surface temperature $T_{\mathrm{s}}$ ($\mathrm{K}$):$${\mathit{\rho }}_{\mathrm{s}}=\mathrm{149.2}+\mathrm{6.84}{\mathrm{ws}}_{\mathrm{10}}+\mathrm{0.48}{T}_{\mathrm{s}},$$with minimum–maximum values of 200–400 $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{3}}$. This parameterisation was defined so that the simulated density of the first 50 cm of snow fits observations collected over the AIS (see Fig. S1, with the snow density database detailed in Table S1 in the Supplement).

• b.

  The aerodynamic roughness length $z_{\mathrm{0}}$ is computed as a function of the air temperature, as proposed in . The parameterisation was tuned so that $z_{\mathrm{0}}$ fit the observed seasonal variation between high ($>\mathrm{1}$ $\mathrm{mm}$) summer and lower (0.1 $\mathrm{mm}$) winter values in coastal Adélie Land, for air temperatures above $-\mathrm{20}$ ${}^{\circ }$$\mathrm{C}$. For lower temperatures, $z_{\mathrm{0}}$ is kept constant and set to 0.2 $\mathrm{mm}$, in agreement with observed $z_{\mathrm{0}}$ values on the Antarctic plateau e.g.;

• c.

  As in , the MAR drifting snow scheme is not activated because this scheme was sensitive to parameter choices . An updated version of the drifting snow scheme is currently being developed and evaluated for application at the scale of the whole ice sheet.

We compare MAR results over the AIS to the latest outputs of the regional atmospheric model RACMO2 version 2.3p2 , called RACMO2 hereafter, using a horizontal resolution of 27 $\mathrm{km}$, a vertical resolution of 40 atmospheric levels, and a topography based on the digital elevation model from . This regional model is developed by the Royal Netherlands Meteorological Institute (KNMI) and has subsequently been adapted for modelling the Antarctic climate and its SMB . It includes a drifting snow scheme , an albedo routine with prognostic snow grain size , and a multilayer snow model computing melt, percolation, refreezing, and run-off .

MAR and RACMO2 models were developed independently. We will not detail here the many physical parameterisation differences between both RCMs, but we will later highlight some of them we show having a significant impact on the modelled SMB.

Forcing reanalyses

Regional atmospheric models are forced by atmospheric fields at their lateral boundaries (pressure, wind, temperature, humidity), at the top of the troposphere (temperature, wind), as well as by sea surface conditions (sea ice concentration, sea surface temperature) every 6 h. Consequently, regional atmospheric models add details and physics to the forcing model in the middle and lower troposphere and at the land or iced surface, whereas large-scale circulation patterns are driven by the forcing fields. We forced MAR with three reanalyses over Antarctica in order to evaluate the uncertainty in the simulated surface climate arising from the uncertainty in the assimilation systems: the European Centre for Medium-Range Weather Forecasts “Interim” re-analysis hereafter ERA-Interim, resolution $\sim \mathrm{0.75}$${}^{\circ }$, i.e. $\sim \mathrm{50}$ $\mathrm{km}$ at 70${}^{\circ }$ S;, the Modern-Era Retrospective Analysis for Research and Applications version 2 hereafter MERRA-2, resolution $\sim \mathrm{0.5}$${}^{\circ }$;, and the Japanese 55-year Reanalysis from the Japan Meteorological Agency hereafter JRA-55, resolution $\sim \mathrm{1.25}$${}^{\circ }$;.

The regional atmospheric model RACMO2 is forced by ERA-Interim. We focus our study to the period 1979–2015, as reanalyses are known to be unreliable before 1979, when satellite sounding data started to be assimilated .

Observations

SMB observations and sectors of strong SMB gradients

We use SMB observations of the GLACIOCLIM-SAMBA dataset detailed in and updated by . This dataset is an update of the one assembled by following the quality-control methodology defined by . It includes 3043 reliable SMB values averaged over more than 3 years. We add accumulation estimates from , retrieved over the Amundsen Sea coast (Marie Byrd Land) with an airborne-radar method combined with ice-core glaciochemical analysis.

The first-order feature of the Antarctic SMB is a strong coastal–inland gradient, with mean values ranging from typically greater than 500 $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$ at the ice sheet margins to about 30 $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$ in the dry interior plateau Fig. a; see also, e.g.. As observations only cover 5 % of MAR grid cells over the ice sheet, we divide that sparse observation dataset into 10 sectors detailed in Table  and shown in Fig. . Six of them are stake transects with a stake every $\sim \mathrm{1.5}$ $\mathrm{km}$, which have been proven very valuable for evaluating modelled SMB . The four other sectors are composed of more scattered observations covering large elevation ranges (Victoria Land, Dronning Maud Land, and Ross Ice Shelf–Marie Byrd Land).

Sectors extracted from the GLACIOCLIM-SAMBA database.

Sector name	Sector type	No. of	No. of	Year range	Elevation range ($\mathrm{m}$ a.s.l.)	Ref.
		obs.	grid cells
Marie Byrd Land	Radar transects	6615	57	1980–2009	973–1873	${}^{\mathrm{1}}$
Ross–Marie Byrd Land	Scattered	72	51	1950–1991	37–1995	${}^{\mathrm{2},\mathrm{3},\mathrm{4}}$
Victoria Land	Scattered	60	40	1951–2006	1804–3240	${}^{\mathrm{5},\mathrm{6},\mathrm{7}}$
Dumont d'Urville–Dome C	Transect	116	24	1955–2010	633–3240	${}^{\mathrm{5},\mathrm{8},\mathrm{9},\mathrm{10}}$
Law Dome–Wilkes Land	Transect	382	32	1973–1986	801–2232	${}^{\mathrm{11}}$
Zhongshan–Dome A	Transect	583	40	1994–2011	1031–4081	${}^{\mathrm{12},\mathrm{13}}$
Mawson–Lambert Glacier	Transect	515	36	1990–1995	1883–2924	${}^{\mathrm{14}}$
Syowa–Dome F	Transect	507	38	1955–2010	584–3803	${}^{\mathrm{15}}$
Princess Elisabeth	Transect	58	6	2009–2012	47–1071	${}^{\mathrm{16}}$
Dronning Maud Land	Scattered	376	104	1955–2008	1753–3741	${}^{\mathrm{17},\mathrm{18},\mathrm{19},\mathrm{20}}$

${}^{\mathrm{1}}$ . ${}^{\mathrm{2}}$ . ${}^{\mathrm{3}}$ . ${}^{\mathrm{4}}$ . ${}^{\mathrm{5}}$ . ${}^{\mathrm{6}}$ . ${}^{\mathrm{7}}$ . ${}^{\mathrm{8}}$ . ${}^{\mathrm{9}}$ . ${}^{\mathrm{10}}$ . ${}^{\mathrm{11}}$ . ${}^{\mathrm{12}}$ . ${}^{\mathrm{13}}$ . ${}^{\mathrm{14}}$ . ${}^{\mathrm{15}}$ . ${}^{\mathrm{16}}$ GLACIOCLIM-BELARE. ${}^{\mathrm{17}}$ . ${}^{\mathrm{18}}$ . ${}^{\mathrm{19}}$ . ${}^{\mathrm{20}}$ .

Model–observation comparison method

RACMO2 outputs are bilinearly interpolated to the $\mathrm{35}\times \mathrm{35}$ $\mathrm{km}$ MAR grid. For each SMB observation, we consider the four surrounding MAR grid cells, from which we eliminate ocean grid cells. We also eliminate surrounding grid cells with an elevation difference with the observation greater than 200 $\mathrm{m}$ (missing elevation of observation is set to Bedmap2 elevation at 1 $\mathrm{km}$ resolution). Finally, we bilinearly interpolate model values of the remaining grid cells at the observation location (see schematic in Fig. S2).

As we restrict our modelling study to the 1979–2015 period, we only consider observations beginning after 1950. For observations beginning after 1979, we time-average model outputs for the same period as the observation. We keep observations beginning before 1979 only if they cover more than 8 years, and in this case we compare the observed value with the modelled value time-averaged for 1979–2015.

In a last step, we average-out the kilometre-scale variability in the observed SMB by binning point values onto grid cells. For each grid cell containing multiple observations, we average all observations contained in the grid cell weighted by the time span of observations, and in the same way we weight-average the modelled values interpolated to observation locations. This way, we obtain consistent observed and modelled averaged values on grid cells.

We discard 66 observations beginning before 1979 and spanning less than 8 years. We also discard 12 observations for which the four surrounding grid cells fall in ocean and seven observations located at specific topographic features for which none of the four surrounding grid cell has an elevation difference of less than 200 $\mathrm{m}$ with respect to the actual location. After this, we retain 559 model–observation comparisons.

MAR SMB for the period 1979–2015: (a) mean annual SMB, with coloured dots showing the observed SMB values (shared colour scale); (b) standard deviation of annual SMB; (c) standard deviation divided by mean annual SMB; (d) difference between MAR and observed SMB on MAR grid cells, following the methodology detailed in Sect. . Magenta dots in panels (b) and (c) show the location of SMB observations. Solid grey lines are contours of surface height every 1000 $\mathrm{m}$ a.s.l. Latitude circles are $-\mathrm{60}$, $-\mathrm{70}$, and $-\mathrm{80}$${}^{\circ }$ S, and longitude lines are from 145${}^{\circ }$ W to 145${}^{\circ }$ E by steps of 45${}^{\circ }$.

[Figure omitted. See PDF]

Modelled vs. observed SMB for sectors and transects as detailed in Table . RACMO2 outputs are bilinearly interpolated to the MAR grid. SMB values are first averaged on MAR grid cells (Sect. ) then along a chosen grid direction (Fig. S2) or by elevation bins. Distance along the transect starts at the coast. Uncertainty of observed SMB (grey shaded area) is the standard deviation of observations contained in each grid cell (sub-grid variability), estimated as a function of the mean observed SMB (see Fig. S3). Despite SMB values corresponding to grid cell averages, we display one marker for each observation, with the $x$ axis corresponding to the observation location along the transect or elevation. For observed SMB plots, markers with white faces are for bins containing fewer than 10 observations and black faces for bins containing more than 10 observations. Magenta bands mark grid cells in which more than 15 % of precipitation sublimates in the katabatic layers according to . The map shows the main Antarctic basins: Antarctic Peninsula in purple, West Antarctic ice sheet in green, and East Antarctic ice sheet in orange. Ice shelves are mapped in blue and grounded islands in red, and the blue line shows the location of the grounding line.

[Figure omitted. See PDF]

Results

Evaluation of the modelled SMB

The large spatial Antarctic SMB gradients, shown in Fig. a as modelled by MAR forced by ERA-Interim for the period 1979–2015, coincide with a strong interannual variability (Fig. b), expressed by a standard deviation of $\sim \mathrm{22}$ % of the mean SMB on average over the ice sheet (Fig. c). MAR SMB shows no systematic spatial bias (Fig. d), with a mean bias of 6 $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$ (4 % of the mean observed SMB), as well as a very strong correlation with the observed SMB ($R^{\mathrm{2}}=\mathrm{0.83}$, $p$ value $<\mathrm{0.01}$, computed on the logarithm of SMB values, as SMB distributions are log-normal). RACMO2 shows similar performance (mean bias of $-\mathrm{3}$ $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$, $R^{\mathrm{2}}=\mathrm{0.86}$, computed on the logarithm of SMB as well).

The model–observation comparison by sectors (Fig. ) reveals a good representation of the coast-to-plateau SMB gradients by both RCMs. MAR and RACMO2 are in good agreement despite MAR not including drifting snow processes whereas RACMO2 does, except in Ross–Marie Byrd Land and in Victoria Land where MAR simulates larger SMB than RACMO2. Another noticeable result is that MAR forced by ERA-Interim, JRA-55, and MERRA-2 gives very similar results for the SMB spatial pattern, not only at the observation locations (Fig. ) but also at the ice sheet scale (comparisons of MAR SMB for different forcing reanalyses are shown in Fig. S4, with colour map scales 10 times smaller than in Fig. S5 in which MAR is compared to RACMO2). This is why we focus on MAR forced by ERA-Interim in the following.

We find no significant differences in the SMB simulated by MAR and RACMO2 when integrated over the ice sheet or its major basins (Table ). SMB is driven by snowfall amounts, which are more than 10 times larger than other SMB components. Snow sublimation in RACMO2 is the sum of sublimation at the surface of the snowpack and of drifting snow sublimation and is approximately 50 % larger than in MAR, which only includes surface snow sublimation. However, surface snow sublimation alone is almost 2 times larger in MAR than in RACMO2 (Table  and spatial patterns shown in Fig. S6), which we investigate in the next section. Modelled surface melt is less than half of the sublimation amount; however liquid water almost entirely refreezes into the snowpack in both models (maps of MAR- and RACMO2-modelled melt amounts are shown in Fig. S7). Temporal variability in the SMB and its components is fully driven in both RCMs by the forcing reanalyses and are therefore strongly correlated with each other (time series shown in Fig. S8). We do not elaborate on the SMB temporal variability here as this aspect will be further detailed in a forthcoming study.

Antarctic integrated SMB on average for 1979–2015 $\pm$ 1 standard deviation of annual values, in gigatonnes per year. Antarctic ice sheet (AIS) and basin geometry are based on Rignot basins , shown in Fig. . RACMO2 is bilinearly interpolated on the MAR grid and the same mask is applied to both models, with area given for this mask. SMB is computed as follows: MAR SMB $=$ snowfall $+$ rainfall $-$ surface snow sublimation $-$ run-off; RACMO2 SMB $=$ snowfall $+$ rainfall $-$ surface snow sublimation $-$ drifting snow sublimation $-$ drifting snow transport $-$ run-off.

Basin	Area (10${}^{\mathrm{6}}$ km${}^{\mathrm{2}}$)	Component ($\mathrm{Gt}$ ${\mathrm{yr}}^{-\mathrm{1}}$)	MAR (ERA-Interim)	RACMO2 (ERA-Interim)
Total AIS	$\mathrm{13.41}$	SMB	$\mathrm{2200}\pm \mathrm{115}$	$\mathrm{2177}\pm \mathrm{107}$
w/o peninsula		Snowfall	$\mathrm{2306}\pm \mathrm{111}$	$\mathrm{2339}\pm \mathrm{107}$
		Rainfall	$\mathrm{6}\pm \mathrm{1}$	$\mathrm{2}\pm \mathrm{1}$
		Surface snow sublimation	$\mathrm{111}\pm \mathrm{10}$	$\mathrm{57}\pm \mathrm{4}$
		Drifting snow sublimation	–	$\mathrm{101}\pm \mathrm{5}$
		Drifting snow transport	–	$\mathrm{5}\pm \mathrm{0}$
		Run-off	$\mathrm{1}\pm \mathrm{1}$	$\mathrm{1}\pm \mathrm{1}$
		Melt	$\mathrm{40}\pm \mathrm{20}$	$\mathrm{68}\pm \mathrm{30}$
Total AIS	$\mathrm{13.83}$	SMB	$\mathrm{2517}\pm \mathrm{111}$	$\mathrm{2516}\pm \mathrm{105}$
Grounded AIS	$\mathrm{12.04}$	SMB	$\mathrm{1923}\pm \mathrm{100}$	$\mathrm{1857}\pm \mathrm{94}$
w/o peninsula		Snowfall	$\mathrm{1995}\pm \mathrm{97}$	$\mathrm{1987}\pm \mathrm{94}$
		Surface snow sublimation	$\mathrm{77}\pm \mathrm{8}$	$\mathrm{39}\pm \mathrm{3}$
		Drifting snow sublimation	–	$\mathrm{87}\pm \mathrm{4}$
Grounded AIS	$\mathrm{12.27}$	SMB	$\mathrm{2120}\pm \mathrm{99}$	$\mathrm{2068}\pm \mathrm{93}$
Grounded East AIS	$\mathrm{9.77}$	SMB	$\mathrm{1170}\pm \mathrm{89}$	$\mathrm{1121}\pm \mathrm{80}$
		Snowfall	$\mathrm{1245}\pm \mathrm{87}$	$\mathrm{1225}\pm \mathrm{82}$
		Surface snow sublimation	$\mathrm{77}\pm \mathrm{6}$	$\mathrm{34}\pm \mathrm{3}$
		Drifting snow sublimation	–	$\mathrm{66}\pm \mathrm{4}$
Grounded West AIS	$\mathrm{2.11}$	SMB	$\mathrm{675}\pm \mathrm{62}$	$\mathrm{643}\pm \mathrm{62}$
		Snowfall	$\mathrm{675}\pm \mathrm{61}$	$\mathrm{668}\pm \mathrm{62}$
		Surface snow sublimation	$\mathrm{1}\pm \mathrm{3}$	$\mathrm{4}\pm \mathrm{1}$
		Drifting snow sublimation	–	$\mathrm{20}\pm \mathrm{2}$
Grounded islands	$\mathrm{0.16}$	SMB	$\mathrm{78}\pm \mathrm{7}$	$\mathrm{93}\pm \mathrm{8}$
Grounded peninsula	$\mathrm{0.23}$	SMB	$\mathrm{198}\pm \mathrm{26}$	$\mathrm{211}\pm \mathrm{27}$

Drifting snow transport features

Annual mean 10 m wind speed, curvature of elevation, and modelled SMB minus observed SMB for the four long transects: (a) Low Dome–Wilkes Land, (b) Zhongshan–Dome A, (c) Mawson–Lambert Glacier, and (d) Syowa–Dome F. Blue lines and colour shading are for MAR (ERA-Interim) outputs and red lines are for RACMO2 (ERA-Interim) outputs. Values are computed as in Fig. 2.

[Figure omitted. See PDF]

(a) Difference in SMB by grid cell ($\mathrm{\Delta }$SMB) between MAR (ERA-Interim) and observations for four transects (Law Dome–Wilkes Land, Zhongshan–Dome A, Mawson–Lambert Glacier, and Syowa–Dome F) vs. surface curvature on MAR grid. Curvature is shifted by $\pm$ 1 grid cell according to the maximum correlation with $\mathrm{\Delta }$SMB (Fig. S10). Linear regression through the origin is plotted with a pink dashed line. We excluded the regression of two outliers (dots with black outlines) and seven data for which MAR annual 10 $\mathrm{m}$ wind speed is lower than 7 $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$ (squares with black outlines). (b) Estimate of mean annual drifting snow transport based on a scaling of the curvature: drifting snow transport ($\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$) $=\mathit{\alpha }$ (${\mathrm{10}}^{\mathrm{6}}$ $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{1}}$ ${\mathrm{yr}}^{-\mathrm{1}}$) $\times$ curvature (${\mathrm{10}}^{-\mathrm{6}}$ ${\mathrm{m}}^{-\mathrm{1}}$), with $\mathit{\alpha }=\mathrm{0}$ $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{1}}$ ${\mathrm{yr}}^{-\mathrm{1}}$ for wind speed lower than 5 $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$, $\mathit{\alpha }=\mathrm{3700}$ ${\mathrm{10}}^{\mathrm{6}}$ $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{1}}$ ${\mathrm{yr}}^{-\mathrm{1}}$ for wind speed greater than 9 $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$, and $\mathit{\alpha }$ linearly increasing as a function of wind speed in between, around the 7 $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$ wind speed threshold. Wind speed is the annual mean of 10 m wind speed modelled by MAR (ERA-Interim). Coloured dots show the difference between MAR SMB and observed SMB with the same colour scale. (c) Mean annual drifting snow transport flux in RACMO2 on average for 1979–2015 ($\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$). Coloured dots show the difference between MAR SMB and observed SMB with the same colour scale.

[Figure omitted. See PDF]

Local fluctuations of the observed SMB around the smooth modelled SMB gradients are apparent along the four stake transects covering more than 500 $\mathrm{km}$: Law Dome–Wilkes Land, Zhongshan–Dome A, Mawson–Lambert Glacier, and Syowa–Dome F. We related these fluctuations to drifting snow transport. Indeed, the snow eroded from the snowpack is loaded into the atmosphere, where it can sublimate and be transported by the wind. Katabatic winds blowing on the surface of the ice sheet result from the downslope gravity flow of cold, dense air. As a consequence, the surface wind divergence, which drives the snowdrift mass transport, is strongly related to the curvature of the topography, and both have similar spatial patterns (shown in Fig. S9). This is because slopes becoming steeper (crests, positive curvature) will lead to wind speed acceleration (positive wind divergence), and thus to drifting snow export (mass loss), whereas slopes becoming more gentle (valleys, negative curvature) will lead to wind speed deceleration (negative wind divergence), and thus to drifting snow deposit (mass gain).

To test our hypothesis, we computed the mean curvature of the MAR $\mathrm{35}\mathrm{km}\times \mathrm{35}$ $\mathrm{km}$ elevation field. In Fig. , we notice that both RCMs commonly exhibit an excess of accumulation on crests and a deficit of accumulation in valleys, in the range of $\pm \mathrm{40}$ $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$. To quantify this curvature effect, we correlate MAR SMB bias ($\mathrm{\Delta }$SMB) with the curvature. For each transect, we apply a constant shift of $\pm$ one grid cell to the curvature in order to find the maximum correlation with $\mathrm{\Delta }$SMB. For three out of the four transects, we find only one shift for which the correlation is significant, and for the remaining transect (Syowa–Dome F) we find no significant correlation (Fig. S10). The sign and the amplitude of those shifts are in line with curvature being used as a proxy for wind divergence, as they are consistent with the Coriolis wind deflection westward of the topography gradient (detailed in Fig. S11). After applying those shifts, we find that the difference between modelled and observed SMB ($\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$) is scaled to approximately $\mathrm{3700}\pm \mathrm{1100}$ (${\mathrm{10}}^{\mathrm{6}}$ $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{1}}$ ${\mathrm{yr}}^{-\mathrm{1}}$) times the curvature (${\mathrm{10}}^{-\mathrm{6}}$ ${\mathrm{m}}^{-\mathrm{1}}$), with a significant relationship ($R^{\mathrm{2}}=\mathrm{0.41}$, Fig. a), when the mean annual 10 $\mathrm{m}$ wind speed (ws${}_{\mathrm{10}}$) is greater than 7 $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$. For lower wind speed (ws${}_{\mathrm{10}}<\mathrm{7}$ $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$), we no longer observe any relationship between model bias in SMB and curvature (horizontally aligned squares in Fig. a). This is consistent with the drifting snow transport process, which requires the wind speed to reach threshold values for the erosion to be initiated .

Hence, a large part of the discrepancies between modelled and observed SMB are explained by surface curvature when wind speed is sufficiently high, which we relate to the unresolved drifting snow transport in MAR. We are able to catch the drifting snow transport signal because drifting snow sublimation is negligible for the four studied transects, as they are located at high elevation, above 2000 $\mathrm{m}$ above sea level (a.s.l.), where the cold atmosphere has a low capacity to be loaded with moisture (see detailed analysis in Fig. S12). The moisture holding capacity of the atmospheric boundary layer is an upper bound for drifting snow sublimation and quickly tends to zero when the mean air temperature decreases below $-\mathrm{30}$${}^{\circ }$C, which is the case along most of the transects, whereas the amplitude of observed SMB fluctuations around the smooth SMB gradient is independent of the temperature (Fig. S13).

Consequently, we propose that drifting snow transport fluxes (ds${}_{\mathrm{tr}}$) not resolved by MAR can be estimated as a scaling of curvature depending on wind speed: ds${}_{\mathrm{tr}}=\mathit{\alpha }\left({\mathrm{ws}}_{\mathrm{10}}\right)\cdot \mathrm{curvature}$ (Fig. b). The scaling factor $\mathit{\alpha }\left({\mathrm{ws}}_{\mathrm{10}}\right)$ depends on wind thresholds to simulate the transition between no drifting snow transport for low wind speed ($\mathit{\alpha }=\mathrm{0}$ for ws${}_{\mathrm{10}}<\mathrm{5}$ $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$) and drifting snow transport scaled to curvature for high wind speed ($\mathit{\alpha }=\mathrm{3700}$ ${\mathrm{10}}^{\mathrm{6}}$ $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{1}}$ ${\mathrm{yr}}^{-\mathrm{1}}$ for $w{s}_{\mathrm{10}}>\mathrm{9}$ $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$), with a linearly increasing scaling factor between $\mathrm{5}$ and $\mathrm{9}$ $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$ for a smooth transition around the $\mathrm{7}$ $\mathrm{m}$ ${\mathrm{s}}^{-\mathrm{1}}$ wind threshold defined above. That estimate of drifting snow transport fluxes shows little sensitivity to the choice of the wind thresholds and of the scaling factor (see fluxes summed over the ice sheet for different thresholds and scaling factors in Table S2). The spatial pattern of drifting snow transport we obtain is comparable to the one simulated by RACMO2 (Fig. c), except that it gives fluxes more than 3 times larger than in RACMO2 (Table S2, and note the different colour map scales between Fig. b and c). The drifting snow transport estimate consists in a redistribution of mass with negligible net mass loss over the AIS (total AIS mass gain of $\sim \mathrm{75}$ $\mathrm{Gt}$ ${\mathrm{yr}}^{-\mathrm{1}}$ and total AIS mass loss of $\sim \mathrm{80}$ $\mathrm{Gt}$ ${\mathrm{yr}}^{-\mathrm{1}}$; see Table S2).

Our drifting snow transport estimate gives a good constraint for drifting snow fluxes above 2000 $\mathrm{m}$ a.s.l., where low temperatures induce negligible atmospheric sublimation. As drifting snow transport is proportional to the amount of snow in suspension in the atmosphere, quantifying this flux also enables us to constrain the amount of snow eroded from the snowpack to the atmosphere, which drives drifting snow sublimation fluxes at lower elevation. This is of importance as drifting snow sublimation is a much larger mass sink than drifting snow transport over the whole ice sheet but is still poorly constrained because observations are very scarce below 2000 $\mathrm{m}$ a.s.l. where it occurs.

Drifting snow sublimation included in RACMO2 and not in MAR moistens the surface atmospheric layers, consequently reducing the sublimation at the surface of the snowpack. This might explain the stronger surface snow sublimation in MAR than in RACMO2 (Table  and Fig. S6). However, drifting snow sublimation is a potentially larger mass sink than surface snow sublimation, as drifting snow particles are continuously ventilated and fully exposed to the ambient air. Consequently, by accounting for drifting snow in MAR we expect that the drifting snow sublimation mass sink could be enhanced at the expense of surface snow sublimation at the ice sheet margins.

Sublimation of precipitation in low-level atmosphere

The four maps show mass fluxes in $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$ for the year 2015. (a) Difference in SMB between MAR and RACMO2. Blue lines delimitate areas where the SMB difference is 30 % greater than MAR SMB, with solid lines when MAR is greater than RACMO2 and dashed lines when MAR is lower than RACMO2. (b) Same as (a) but for the snowfall amounts at the ground. (c) Same as (a) but for the sublimation of precipitation in the atmospheric layers. (d) Same as (a) but for the maximum snowfall amount (equal to ground snowfall plus atmospheric sublimation). Locations of transects A1–A2 and B1–B2 extracted in Fig.  are shown in panels (b) and (d).

[Figure omitted. See PDF]

As described above, MAR and RACMO2 RCMs forced with ERA-Interim simulate similar spatial patterns for SMB compared to observations (Fig. ). However, at the ice sheet scale, MAR and RACMO2 SMB show regional discrepancies (shown in Fig. a for 2015, and similar to the 1979–2015 mean shown in Fig. S5a), which are primarily the result of differences in simulated snowfall rates (Figs. b and S5b). We notice that areas where MAR snowfall is much lower than RACMO2 snowfall (Fig. b, blue dashed lines) coincide almost exactly with the pattern of precipitation that is able to sublimate in the low-level atmosphere according to . In that study, the amount of atmospheric sublimation is quantified for the year 2015 using atmospheric modelling constrained with precipitation radar observations. Atmospheric sublimation happens because the katabatic surface air flux, moving from highly elevated inland plateau toward sea level, is subject to adiabatic compression when it moves downslope. This compression induces an increase in air temperature, which reduces relative humidity and drives sublimation rates in the lower troposphere ($\sim$ first 1000 m above the ground), enhanced in the katabatic channels at the ice sheet margins.

To deepen this analysis, we re-ran MAR for the year 2015 in order to save the full atmosphere snowfall fields. From the daily 3-D snowfall amounts, we derived the atmospheric sublimation amount from the difference between the maximum snowfall and the ground snowfall in each atmospheric column, as in . The same was done for RACMO2. We find that the atmospheric sublimation simulated by MAR (363 Gt for the year 2015 over the grounded ice sheet) is higher than estimated in (299 Gt after interpolation on the same mask) and much higher than simulated by RACMO2 (128 Gt, Fig. c). A major difference between MAR and RACMO2 is the advection of precipitation in the atmosphere: in MAR, precipitating particles are explicitly advected through the atmospheric layers until they reach the surface, while in RACMO2, precipitation is added to the surface without horizontal advection, and is able to interact with the atmosphere only in a single time step (6 min in this simulation). Consequently, atmospheric sublimation is likely to be underestimated in RACMO2.

We conclude, in agreement with , that atmospheric sublimation is a major mass sink at the ice sheet margins in MAR, as for the year 2015 it represents 16 % of the snowfall loaded on the grounded ice sheet 12 % in and 26 % for areas below 1000 $\mathrm{m}$ a.s.l. 17 % in.

It is noticeable that very few SMB observations are available in areas where identify low-level sublimation, marked by magenta bands in Fig. . Except for Ross–Marie Byrd Land, the only other areas where low-level sublimation is greater than 15 % of the total precipitation as defined by are close to Dumont d'Urville (coastal Adélie Land) and to Syowa (coastal Dronning Maud Land). In those areas the SMB amount is indeed larger in RACMO2 than in MAR and in observations. Both RCMs overestimate SMB around 2000 $\mathrm{m}$ a.s.l. in Dronning Maud Land and in Ross–Marie Byrd Land (Fig. ), which could indicate katabatic channels not being resolved enough by the topography of the models.

Precipitation formation and advection

MAR- and RACMO2-simulated fields for the year 2015, extracted with a bilinear interpolation for (left) transect A1–A1 and (right) transect B1–B2 (locations shown in Fig. b and d). Each panel shows MAR fields (blue lines) and RACMO2 fields (red lines) for (a) surface height, in metres above sea level; (b) maximum snowfall amounts, equal to ground snowfall plus atmospheric sublimation, in $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$; and (c) snowfall amounts at the ground, in $\mathrm{kg}$ ${\mathrm{m}}^{-\mathrm{2}}$ ${\mathrm{yr}}^{-\mathrm{1}}$. In (b) and (c), the thick black line is for the difference in snowfall between MAR and RACMO2 (MAR-RACMO2), with green-filled areas when MAR snowfall is larger than RACMO2 snowfall, and brown-filled areas when MAR snowfall is lower than RACMO2 snowfall (same convention as in Fig. ); the dotted lines are for the atmospheric sublimation modelled by MAR (blue) and by RACMO2 (red), negative when it induces a decrease in precipitation; light coloured bands show crests (light blue, curvature of MAR topography greater than $\mathrm{0.005}{\mathrm{10}}^{-\mathrm{6}}$ ${\mathrm{m}}^{-\mathrm{1}}$) and valleys (light yellow, curvature of MAR topography lower than $-\mathrm{0.005}{\mathrm{10}}^{-\mathrm{6}}$ ${\mathrm{m}}^{-\mathrm{1}}$). The thick black arrows show the main 800 $\mathrm{hPa}$ wind direction during cyclonic activity.

[Figure omitted. See PDF]

Differences between MAR and RACMO2 snowfall fields are strongly reduced when considering the maximum snowfall amounts (before sublimation in the low-level atmosphere) rather than the ground snowfall amounts (Fig. b and d). However, MAR snowfall rates generally exceed those simulated by RACMO2, by more than 30 % on the lee side of the West AIS (Marie Byrd Land toward Ross ice shelf), on the lee side of the Transantarctic Mountains (Victoria Land), and close to crests at the ice sheet margins. MAR maximum snowfall rates are lower than simulated by RACMO2 windward of topographic barriers and in valleys at the ice sheet margins. This spatial pattern looks similar to the one obtained in RACMO2 when delaying the conversion of cloud ice–water into snow–rain Fig. 3a of. This change led to both ice and water clouds lasting longer in the atmosphere before precipitating and therefore being advected further towards the ice sheet interior .

For a more in-depth analysis, we extract MAR and RACMO2 snowfall rates on two transects at the ice sheet margins (Fig. ), following the main wind direction during cyclonic activities (locations shown in Fig. b and d). On these transects the observed difference in maximum snowfall between MAR and RACMO2 is largely explained by a phase difference in the snowfall peaks windward of the topographic barriers, with MAR peaking closer to the crests than RACMO2 (Fig. b). This induces a wave-like pattern of precipitation difference strongly related to the shape of the topography, with larger snowfall amounts in MAR than in RACMO2 just windward of crests and lower snowfall amounts in MAR than in RACMO2 around windward valleys. At the ground, lower snowfall in MAR than in RACMO2 in valleys is amplified by low-level atmospheric sublimation, which peaks in katabatic channels (Fig. c).

Observations do not enable us to definitively discriminate one model against the other, but we observe a general tendency for MAR to overestimate accumulation on Ross–Marie Byrd Land and close to ice sheet summits (Dome C, Dome A, Dome F; see Figs. d and ). Close to summits the wind is low, so a missing drifting snow transport process is an unlikely explanation for a positive bias in SMB modelled by MAR (Fig. b). Over the Greenland ice sheet, MAR tends to overestimate ice cores based on accumulation inland while RACMO2 underestimates it .

We conclude that the differences in MAR and RACMO2 snowfall patterns are very likely related to differences in the advection of precipitation inland, which may arise from (i) the different advection of precipitating particles to the ground described in Sect. , (ii) different timing of precipitation formation (cloud–precipitation conversion thresholds), and/or (iii) different dynamical response to the topographic forcing, caused by either different dynamical cores or the different resolutions (the 27 $\mathrm{km}$ resolution in RACMO2 better resolves the ice sheet topography than the 35 $\mathrm{km}$ resolution in MAR).

Discussion and conclusion

In our study, we evaluate new estimates of the Antarctic SMB obtained with the polar RCM MAR ran for the first time for decades-long simulations at the scale of the whole AIS. We use model settings comparable to previous MAR simulations over Greenland but with a specific upper atmosphere relaxation and new surface snow density and roughness length parameterisations. We present simulations of MAR forced by ERA-Interim, JRA-55, and MERRA-2 for the satellite era (1979–2015) in which we can rely on reanalyses products. Remarkably, MAR forced by those three reanalyses gives similar spatial and temporal SMB patterns. We also compare MAR with the latest simulations of the RCM RACMO2 forced by ERA-Interim . We find no significant differences between MAR and RACMO2 SMB when integrated on the AIS and its major basins (Table ).

As the dominant feature of the Antarctic SMB is its strong coast-to-plateau gradient, we extract stake transects and sectors with large elevation ranges from the GLACIOCLIM-SAMBA SMB observational dataset. We show that both RCMs show similar performances when compared to observations, with a good representation of the SMB gradient (Fig. ). But more importantly, we outline and quantify missing or underestimated processes in both RCMs.

Along stake transects, we relate 100 $\mathrm{km}$ scale fluctuations of observations around the smooth modelled SMB pattern to the shape of the ice sheet captured on the $\mathrm{35}\mathrm{km}\times \mathrm{35}$ $\mathrm{km}$ MAR grid. Both RCMs accumulate too much snow on crests, and not enough snow in valleys, as a result of drifting snow transport fluxes not included in MAR and probably underestimated in RACMO2 by a factor of 3 (Fig. ). In the RACMO2.3p2 version used here, the modified drifting snow routine induced almost halved drifting snow transport and sublimation fluxes compared to the previous RACMO2.3p1 version . In a recent study combining satellite observation of drifting snow events and reanalysis products, estimate the drifting snow sublimation to be about $\sim \mathrm{393}$ $\mathrm{Gt}$ ${\mathrm{yr}}^{-\mathrm{1}}$ over the AIS, vs. 181 $\mathrm{Gt}$ ${\mathrm{yr}}^{-\mathrm{1}}$ in RACMO2.3p1 and 102 $\mathrm{Gt}$ ${\mathrm{yr}}^{-\mathrm{1}}$ in RACMO2.3p2 . Consequently, observational constraints from our study and from both tend to confirm that drifting snow transport and sublimation fluxes are likely much larger than previous model-based estimates and need to be (better) resolved and constrained in climate models.

We also point out that MAR generally simulates larger SMB and snowfall amounts than RACMO2 inland, particularly on the lee side of the Transantarctic Mountains and on crests at the ice sheet margins, whereas MAR simulates lower snowfall than RACMO2 windward of mountain ranges and promontories. Sublimation of precipitating particles in low-level atmospheric layers is largely responsible for the significantly lower snowfall rates in MAR than in RACMO2 in valleys at the ice sheet margins. As precipitating snow particles have larger time residence in the atmosphere in MAR than in RACMO2 (Sect. ), amounts of precipitation lost by sublimation in katabatic channels are more than twice as much in MAR as in RACMO2. The remaining spatial differences in snowfall between MAR and RACMO2 are attributed to differences in advection of precipitation, with snowfall particles being likely advected too far inland in MAR.

Atmospheric sublimation represents 429 $\mathrm{Gt}$ ${\mathrm{yr}}^{-\mathrm{1}}$ in MAR over the whole AIS (peninsula excluded) for the year 2015, 89 % of which is lost below 2000 $\mathrm{m}$ a.s.l. and 61 % below 1000 $\mathrm{m}$ a.s.l. This might be of importance for the mass balance of glacier drainage basins SMB minus discharge;, as ice streams are typically channel-shaped areas affected by low-level sublimation of precipitation. Consequently, we note the importance of saving precipitation fluxes in models at least 1300 $\mathrm{m}$ above the ground for comparison with CloudSat products, but ideally at all model levels below 1500 $\mathrm{m}$ above the ground to be able to compute sublimation of precipitation in the low-level atmospheric layers. This will become a standard output in forthcoming MAR simulations.

We expect that accounting for drifting snow in MAR will lead to significant improvements in describing the Antarctic SMB and surface climate, as it will enable (1) a quantification of the drifting snow sublimation mass sink, (2) a more realistic representation of relative humidity and temperature in the boundary layer, and (3) an explicit modelling of the drifting snow transport from crests to valleys. Exploring the impact of horizontal and vertical model resolution on drifting snow estimates and on sublimation of precipitation in katabatic channels will also be of importance as those processes are related to the shape of the ice sheet and to the advection of precipitation in the atmosphere. The accuracy of the topography has to be considered as well, as digital elevation models are in constant improvement over the AIS e.g. and should be regularly updated in climate models.

Python scripts developed for this study as well as all required data are available at 10.5281/zenodo.2548847 (). The last version of MAR is freely distributed at http://mar.cnrs.fr/ (last access: 24 January 2019). Monthly MARv3.6.4 outputs from this study are freely available at 10.5281/zenodo.2547637 (), together with the associated MAR source code. The ECMWF reanalysis ERA-Interim 6-hourly outputs were downloaded from 10.5065/D6CR5RD9 (). The MERRA-2 reanalysis 6-hourly outputs were downloaded from https://disc.gsfc.nasa.gov/datasets/ (). The JRA-55 reanalysis 6-hourly outputs were downloaded from 10.5065/D6HH6H41 ().

The supplement related to this article is available online at: https://doi.org/10.5194/tc-13-281-2019-supplement.

CA set-up the MAR model for Antarctica with several adaptations, performed model simulations and analysed model outputs and observations. CA, AO, XF, and VF designed the study. CA, XF, HG, CA, and CK developed the MAR model and contributed to the MAR set-up and output analyses. XF and HG are the main developers of the MAR model. MRvdB, JMvW, JTML, and WJvdB contributed to RACMO2 output analyses. All authors contributed to discussions in writing this paper.

The authors declare that they have no conflict of interests.

Acknowledgements

We thank Kenichi Matsuoka, Massimo Frezzotti, and the anonymous reviewer for their constructive and insightful comments, which led to a much improved paper. Cécile Agosta performed MAR simulations during her Belgian Fund for Scientific Research (F.R.S.-FNRS) research fellowship. Computational resources have been provided by the Consortium des Équipements de Calcul Intensif (CÉCI), funded by the F.R.S.-FNRS under grant no. 2.5020.11. We acknowledge Jacopo Grazioli for sharing the low-level sublimation product and his expertise on this dataset. We acknowledge Yetang Wang for sharing his updated version of the GLACIOCLIM-SAMBA dataset. We acknowledge Christophe Genthon for fruitful discussions and suggestions. The authors acknowledge the support from Agence Nationale de la Recherche Scientifique for the scientific traverses in Antarctica and the associated research on climate and surface mass balance modelling, projects ANR-14-CE01-0001 (ASUMA) and ANR-16-CE01-0011 (EAIIST). Cécile Agosta acknowledges the support from Fondation Albert 2 de Monaco under the project Antarctic-Snow (2018–2020). Edited by: Kenichi Matsuoka Reviewed by: Massimo Frezzotti and one anonymous referee