Snow dune growth increases polar heat fluxes

Falling snow often accumulates in dunes. These bedforms are found on up to 14% of the surface of Earth, and appear occasionally on other planets. They have been associated with increased heat fluxes and rapid sea ice melting (Petrich et al., 2012; Popović et al., 2018). Their formation, however, is poorly understood (Filhol and Sturm, 2015; Kochanski et al., 2019a; Sharma et al., 2019). Here, we use field observations to show that dune growth is controlled by snowfall rate and wind speed. 5 We then use numerical experiments to generate simulated dune topographies under varied wind and snowfall conditions, and use those to quantify conductive and radiative heat fluxes through snow. Our results show that dune growth leads to decreased snow cover, more variable snow depth, and significant increases in surface energy fluxes. We provide quantitative results that will allow modelers to account for the impact of snow bedforms in snow, sea ice, and climate simulations. In addition, this work offers a starting point for process-based studies of one of the most widespread bedforms on Earth. 10


Introduction
Snow dunes are ubiquitous on sea ice (Petrich et al., 2012;Popović et al., 2018), common in Antarctica (Doumani, 1967;Mather, 1962;Watanabe, 1978), and seasonally widespread across Siberia, Alaska, Canada, Scandanavia, Japan and the Rocky Mountains (Bagnold, 1941;Colony et al., 1998;Filhol and Sturm, 2015;Sturm and Liston, 2003;Cornish, 1902;Filhol et al., 2017;Kobayashi, 1980;Kochanski et al., 2018Kochanski et al., , 2019a. Dunes and other bedforms are often encountered by skiiers, snowmo-15 bilers, and polar travelers. These areas collectively cover roughly 14% of the surface of the Earth, implying that snow dunes could be one of the most widespread bedforms. They are found, however, in sparsely populated areas, and have received far less scientific attention than dunes of sand. Existing literature has shown that snow dunes significantly alter the surface energy fluxes. For example, Sturm and Liston (2003) observed that average heat fluxes through Arctic lakes and tundra were twice as high through dune-covered snow as 20 through smooth snow. Several recent studies have shown that dunes control melt on Arctic sea ice by concentrating heat fluxes in narrow inter-dune areas and exposing dark, snow-free ice to the sun (Petrich et al., 2012;Liston et al., 2018;Popović et al., 2018). Additional literature has shown that other snow bedforms, such as sastrugi, increase surface roughness and turbulent heat transfer (Jackson and Carroll, 1978;Leonard, 2009); snow dunes likely have similar effects.
Finally, there is an extensive literature on the growth of dunes made from sand. Dune growth is an instability that occurs 25 when wind stronger than a critical speed blows over a loose, granular surface (Bagnold, 1941;Werner, 1995;Elbelrhiti et al., 2005). The orientation, shapes, and sizes of the dunes depend on wind direction and sediment supply (Andreotti et al., 2002;Gao et al., 2015), with sediment typically coming from a line or point source upwind of the dune field, such as an eroding rock feature. The dune dynamics may be additionally altered by armoring due to non-uniform grain size (Gao, Xin, Clément Narteau, and Olivier Rozier, 2016) or the erosion of sand particles along the dune field (Jerolmack,D. J.,Reitz,M. D.,Martin,30 R. L., 2011). Snow dunes are therefore a special case of dune growth with (1) sediment of lower density and lower average diameter than sand, (2) more fragile and erodible sediment (Comola et al., 2017) than sand, (3) sediment supplied across a wide area from the atmosphere, rather than sediment supplied from a line or point source, and (4) sediment supply that varies over minutes, hours, and days.
In this paper, we focus on the impacts of snow supply and wind speed on dune growth. We show that, unlike sand dunes, 35 which are formed over months or centuries and often reach stable equilibria, snow dunes evolve rapidly during storms and their shapes and sizes change quickly over time.

Observations of snow accumulation
We monitored bedform growth for three winters in the Colorado Front Range using time-lapse photography. Four types of surfaces formed in newly deposited snow; they are shown in Fig. 1. Fig. 1a shows flat snow, the geometry assumed by most 40 climate and snow science models. Fig. 1b shows thin (< 2 cm) patches. These split, merge, and disappear from minute to minute. They typically have narrow upwind ends and one or two downwind points. Fig. 1c shows two barchan dunes. These are more stable than patches. They are distinguished by sharp crests and downwind avalanche faces, and move downwind at speeds of up to 2.5 m/h (Kochanski et al., 2019a). Finally, Fig. 1d shows a field of transverse snow waves. These have intermittent crests, 3-15 m wavelengths (Filhol and Sturm, 2015), and extend tens or hundreds of meters perpendicular to the 45 wind. Any of these surfaces may be further adorned by small (< 2 cm high) marks, like ripples. Additional imagery and criteria for distinguishing the surface types is given in the supplementary material.

Statistical methods for field data
The results below summarize observations from three seasons of fieldwork on Niwot Ridge, at 3528 m elevation in the Front Range. The ridge is a natural wind-blown snow laboratory: it receives some dozen large snowstorms each winter, which are 50 blown over the tundra by unidirectional west-northwesterly winds. We were therefore able to observe snow accumulation in the absence of confounding factors like variable wind directions. The winter wind blows consistently west-northwest and averages 10.5 m/s with gusts up to 29.7 m/s. The full observations are archived in Kochanski (2018b); representative excerpts are viewable at Kochanski (2018a). Climate records are drawn from Morse and Losleben (2019). Our data (Kochanski, 2018b) are by an order of magnitude the most extensive snow bedform observations collected to date.

55
Three experts labeled our time-lapse videos at ten-minute intervals, with labels: 'flat', 'patches', 'dunes', 'waves', 'unidentifiable' and 'not depositional'. Conflicting labels were resolved in favor of the majority. For this study, which focuses on dry, newly deposited snow, we excluded surfaces that had melted (temperature > −1 o C, the most predictive temperature threshold in our data (Kochanski et al., 2018)) and 'not depositional' observations. The resulting data set was 24% flat snow, 24% patches, 21% dunes, 28% waves, and 3% unidentifiable. These examples, excluding unidentifiable surfaces, are shown as functions of 60 wind speed and time in Figs. 1b and c.
We then built softmax classifiers, a method for separating classes of data by logistic regression, to identify the most likely surface type as a function of weather conditions. We tested classifiers that make predictions using combinations 1-3 of weather variables, from: time, wind speed, gust speed (highest wind speed recorded over each 10 minute interval); and time-averaged wind speed, wind stress, and wind power since snowfall. We excluded wind direction, humidity, and temperature, three vari-65 ables that had negligible predictive power in previous work with this data (Kochanski et al., 2018).
We quantified the uncertainty in the weather variables using the Monte-Carlo method described in Kochanski et al. (2018), generating 20 data points from each 10 minute video observation, for a total of 11,202 training examples (not an even multiple of 20 because of uncertainty in which observations exceeded melting temperature) representing 93 hours of observations. The classifiers were trained on all labeled examples and were validated by bootstrapping (Efron and Tibshirani, 1985). The  The presence of flat snow was predicted better by wind speed than by any other variable(s). Dune and wave presence were predicted comparably well by time and wind speed, time and time-averaged wind stress, or time and time-averaged wind power.
We use the wind speed classifier in the rest of this text for simplicity and consistency with simulations.

Fieldwork results
Immediately after snowfall, most snow was flat at wind speeds slower than (4.3 ± 0.9) m/s measured 7.5 m above the surface ( Fig. 2a) and most snow was bedform-covered if wind speeds were higher. We hereafter use 4.3 m/s as the value for the critical wind speed, u c , required to initiate snow saltation and bedform growth at our field site. This value is within the range of previous measurements of u c for dry snow (Li and Pomeroy, 1997). At wind speeds from 4.3 m/s to 12.0 m/s most snow 80 accumulated in patches, and at wind speeds higher than (12.0 ± 2.7) m/s most snow accumulated in dunes and waves.
Over time, the probability of observing flat and dune-covered snow surfaces decreased while the probability of observing patches and snow-waves increased (Fig. 2b). The flat surfaces were routinely remodeled by strengthening winds; every surface we observed experienced winds higher than 4.3 m/s within 1.5 days of snowfall. Once bedforms grew, the surface did not become flat again unless the bedforms were buried by future storms.

85
For the next section, we employ a numerical model to explore how snow supply governs transitions between patches, waves, and dunes. To quantify the 3D dynamics of dune growth in a measurable environment, we simulated accumulation on an initially-flat surface with constant snowfall (s) and constant wind speeds (u). We simulated snow dune growth using Rescal-snow (Kochanski et al., 2019b), an adaptation of sand dune model ReSCAL (Rozier and Narteau, 2014), a community-standard model that has been used in several ground-breaking studies of sand dunes (Gao et al., 2015;Zhang et al., 2012). Rescal-snow is a stochastic 3D cellular automaton that models granular physics through the interactions of neighboring cells, each of which represents a 95 parcel of snow or air. The wind is represented through a lattice gas model. Rescal-snow is fully described by Kochanski et al. (2019b).
We ran fifty simulations. Each modeled 3-7 days of constant snowfall (0.20, 0.63, 2.0, 6.3, or 20 mm/h) blown by constant wind (u/u c = 1, 1.005, 1.05, 1.1, 1.3, 2.0, 2.6, 3.2, 3.9 or 4.5) over an initially flat, bare surface. We modeled a 1200l 0 ×150l 0 domain (about 60m × 7m; unit l 0 is discussed below); this is several times the longest observed bedform wavelength. Snapshots 100 in Fig. 3a, b show a 150l 0 × 190l 0 (about 7.5m × 9.5m) excerpt. Wind was driven and measured along the upper boundary at height 100 l 0 (about 5 m). Grains appeared stochastically along the upper boundary to simulate snowfall, and exited the simulation by blowing through the open downwind boundary. We used periodic wind-parallel boundaries to emulate a wider domain. 105 We calibrated the length (l 0 ), time (t 0 ) and stress (τ 0 ) scales in the non-dimensional model by matching three simulated quantities to their observable counterparts: the maximum unstable wavelength of a granular bed; the wind speed u/u c ; and the saturated flux of wind-blown snow grains. We assume the following conditions, typical of our field site: snow grains of diameter  Fig. 2a. We assume the atmosphere is stable, so that the friction velocity u * = uκ/ ln(z/z 0 )

Model calibration
for wind velocity u measured at height z and Von Kármán constant κ = 0.4. We then match several simulated lengths, times, and stresses with their real-world counterparts in order to calibrate the model. The scaling calculations are presented fully in supplementary § S1.3.
The calibration is the dominant source of uncertainty in this simulation. We are able to match lengths and heights to within 115 ±50 %; times to within ±80 %; and wind speeds/stresses to within ±20%. To minimize the impact of this uncertainty, we present most of our results in non-dimensional forms, e.g. snow cover as a fraction of 100% (f ); roughness length as a fraction of snow depth (σ/h avg ); and heat fluxes through dune-covered snow as a fraction of heat fluxes through uniformly-deep flat snow (q/q flat ). The main exception is the use of results "after 10 cm of snowfall", which to be fully precise should be interpreted as, "after 5-15 cm of snowfall, depending on air density and snow grain size". showed that the variability of the snow and ice albedos and thermal conductivities leads to a 1% variability in the effect of 135 bedforms on conductive heat fluxes, q C /q flat C , and a 40% variability in the effect of bedforms on radiative heat fluxes, q S /q flat S . This result will be unsurprising to snow scientists, as small changes in snow and ice albedo are known to have huge effects on Arctic thermodynamics.
The conductive heat flux q C varies with snow thickness, h(x, y), and the thermal conductivities of snow and ice. For simplicity, we calculate the one-dimensional steady-state heat flux: where A is the sample area and T s − T b is the temperature difference between the surface and ice-ocean interface. Fig. 4c presents q C /q flat C := q C (h s )/q C (h s = h avg ). This quantity is independent of temperature. This calculation neglects lateral heat fluxes, and will therefore err on the side of under-estimating the impact of snow dunes on total heat fluxes.
The shortwave radiation flux q S varies with the snow cover fraction f and the relative albedos of the snow and the surface 145 beneath: where S 0 is the incoming shortwave energy flux. Fig. 4d presents q S /q flat S := q S (f )/q S (f = 1). This quantity is independent of S 0 . Fig. 3a shows how wind speed controls the surface topography after one day of snowfall. Simulated winds at or below the critical speed (u/u c ≤ 1) produced flat snow; gentle winds (1 < u/u c ≤ 2) produced low-lying patches; moderate winds (2 < u/u c ≤ 3.2) mediated a gradual transition from patches to dunes; and high winds (u/u c ≥ 3.2) produced well-defined dunes with slip faces. These results agree quantitatively with our field observations, in which flat surfaces only appeared when 155 u/u c ≤ 1, and dunes became more common than patches when u/u c ≥ 2.5 ± 1.0. The simulations also demonstrated the mechanism by which dunes and waves co-exist at high wind speeds: when simulated dunes touched arms with their neighbors, they briefly formed long waves, then parted. This process resembles the movement of dislocations in a ripple field, played out at a larger length scale. Snowfall rate is fixed at 2 mm/h in a-b. Fig. 3b shows how dunes grew with time given a constant snow supply (2 mm/h) and a constant wind (u/u c = 2). Snow 160 initially (12 h) collected in patches. After 24-36 h, the patches grew into better-defined dunes, then continued to grow (36-72 h) until they overlapped their neighbors and merged into waves. These results largely agree with our field observations, in which patches were the most likely initial bedform at this wind speed, and dunes transitioned into waves after 24-48 h.

Comparison of observed and modeled snow dunes
Figs. 3c-f show how bedform growth altered snow depth distribution and cover fraction. Bedforms grew taller with either stronger winds (Fig. 3c) or more time and total snowfall (Fig. 3d). Dunes produced by high winds were widely separated and 165 covered little of the surface (Fig. 3e). In all cases, the surface cover fraction increased with time (Fig. 3f).
Snow that fell faster (Figs. 3c-f, dashed orange lines) formed shorter dunes that covered the surface more quickly than slow-falling snow (dashed purple lines), and thus covering the surface completely with a smaller volume of fallen snow. This occurred because faster snowfall initiated a larger number of dunes than sparse snowfall, with a correspondingly higher surface area to volume ratio. Bedform size and snow cover fraction both change the surface energy balance, as discussed in § 1. We used our simulated results to quantify the impact of high winds and low snowfall rates on snow depth variability and cover fraction , over a range of conditions observed on Earth (Fig. 4c). We found that bedforms alter the standard deviation of snow depth by up to 80% in the studied range, and the cover fraction by up to 56%. 175 We then used simulated snow topographies to estimate changes in conductive and shortwave energy fluxes through 10 cm of snow atop a 35 cm thick ice floe. We normalized these heat fluxes (Fig. 4d-e) by the heat fluxes through flat snow so that the results are thus independent of temperature and radiation intensity ( § 3.1.3). Fig. 4c shows wind speeds and snowfall rates for seven sites on ice sheets, sea ice, and land, to help readers match the conditions we describe with places they are familiar with. Exact locations, chosen for data availability, are: Savoonga, St.
Snow dune growth increased the conductive heat flow through the ice by up to 34% within the studied range (Fig. 4d). Dune growth had a greater impact on the balance of shortwave radiation, which it increased by up to 94% (Fig. 4e), as the snow dunes left up to 56% of the less-reflective underlying ice exposed. The magnitude of this radiative effect varies by ±34% with reasonable albedo values for snow and ice. The effect would be much greater if the underlying sub-snow were darker than sea 190 ice (e.g. soil, vegetation, most natural surfaces), and moderately smaller if the surface were brighter (e.g. old snow). If this example ice floe were north of the Arctic circle, and had the largest dunes we simulated (grown from a wind speed of 2u c and snow falling at 0.2 mm/h), dune growth would hasten the conductive winter freezing rate of the floe by 34%, and accelerate the radiation-driven melt rate by 94% in spring.

195
Many models of snow cover, sea ice, and the Earth surface include estimates of snow depth, snow cover fraction and surface roughness. The following equations, obtained by statistical fits to our simulation results, allow calculation of those variables accounting for the effects of bedform growth. When the wind speed is small or the snow depth is large, these equations simplify to the forms expected for flat snow.
The snow cover fraction f and standard deviation σ of snow depth obtained from our simulations may be estimated from 200 snowfall rate s (see estimates below), depth of snow fallen h avg , and wind speed u/u c as: where h f = 1700 − 37s, σ max = 140 mm, σ 0 = 110 mm, h σ = 54s + 55, and snowfall rate 0.6 ≤ s ≤ 20 mm/h. All lengths are 205 in millimeters; snowfall rates are in mm/h. The fit approximates our simulation results for f to within ±8% and σ to within ±13 cm.

Conclusion
Flat snow should be expected only where winds are consistently slower than the critical speed for snow transport. Such gentle winds are typical of sheltered forests, but not of tundra, ice sheets, sea ice, or the poles (Fig. 5). Snow dunes have been 210 documented in Japan (Kobayashi, 1980), Alaska (Filhol and Sturm, 2015), Colorado (Kochanski et al., 2019a), Antarctica (Doumani, 1967;Gow, 1965), Greenland (Albert and Hawley, 2002), and on Arctic (Petrich et al., 2012) and Antarctic (Massom et al., 2001) sea ice. Dunes may also grow from snowfall on Mars, Titan, and Pluto. Those three bodies likely have solid precipitation (Head et al., 2005;Perron et al., 2006;Telfer et al., 2018) and certainly have dune-forming winds (Hansen et al., 2011;Lorenz et al., 2006;Telfer et al., 2018). For most snow, on Earth and in the solar system, we need not ask whether dunes