Spatial and temporal development of incipient dunes
Cyril Gadal, Clément Narteau, Ryan C. Ewing, Andrew Gunn, Douglas Jerolmack, Bruno Andreotti, Philippe Claudin
mmanuscript submitted to
Geophysical Research Letters
Spatial and temporal development of incipient dunes
C. Gadal , C. Narteau , R.C. Ewing , A. Gunn , D. Jerolmack , , B.Andreotti , P. Claudin Institut de physique du globe de Paris, Universit´e de Paris, CNRS, Paris, France Department of Geology and Geophysics, Texas A&M University, College Station, USA Department of Earth and Environmental Science, University of Pennsylvania, Philadelphia, USA Department of Mechanical Engineering and Applied Mechanics, University of Pennsylvania, Philadelphia,USA Laboratoire de Physique, ENS - PSL Research University, Universit´e de Paris, CNRS, SorbonneUniversit´e, Paris, France Physique et M´ecanique des Milieux H´et´erog`enes, ESPCI Paris - PSL Research University, Universit´e deParis, CNRS, Sorbonne Universit´e, Paris, France
An edited version of this paper was published by AGU. Copyright 2020 American Geo-physical Union: Gadal, C., Narteau, C., Ewing, R. C., Gunn, A., Jerolmack, D., Andreotti, B.,& Claudin, P. (2020). Spatial and temporal development of incipient dunes. Geophysical Re-search Letters, 47, e2020GL088919. https://doi.org/10.1029/2020GL088919 . Key Points: • Length and time scales of dune formation set the pace of aeolian landscape dy-namics • We provide significant statistics on the temporal and spatial growth of incipientsand dunes extracted from extensive field observations • Data provide a field validation of dune instability theory, and introduce the dis-tance of dune growth as an important length scale
Corresponding author: Cyril Gadal, [email protected] –1– a r X i v : . [ phy s i c s . g e o - ph ] O c t anuscript submitted to Geophysical Research Letters
Abstract
In zones of loose sand, wind-blown sand dunes emerge due the linear instability of a flatsedimentary bed. This instability has been studied in experiments and numerical mod-els but rarely in the field, due to the large time and length scales involved. We exam-ine dune formation at the upwind margin of the White Sands Dune Field in New Mex-ico (USA), using 4 years of lidar topographic data to follow the spatial and temporal de-velopment of incipient dunes. Data quantify dune wavelength, growth rate, and prop-agation velocity and also the characteristic length scale associated with the growth pro-cess. We show that all these measurements are in quantitative agreement with predic-tions from linear stability analysis. This validation makes it possible to use the theoryto reliably interpret dune-pattern characteristics and provide quantitative constraintson associated wind regimes and sediment properties, where direct local measurementsare not available or feasible.
Plain Language Summary
Dunes are the solar system’s ubiquitous landform, arising wherever wind blows overa loose sand bed. An aerodynamic theory for dune formation, which connects grain-scalemovement to emergent dune pattern, has been developed for idealized scenarios. Yet thismodel has never been directly tested in nature, because of the complexities in observ-ing dune formation at the initial stage. Here we report extensive topographic observa-tions of the initiation, growth, and migration of real-world sand dunes. Moreover, we finda surprisingly precise agreement with the idealized aerodynamic theory. This robust con-firmation of the theory for dune formation means that we may estimate wind conditionsin remote places, including other planets, with confidence.
The development of sand dunes, from incipient to mature bedforms, and their evo-lution, through interaction and coarsening processes, involve characteristic time and lengthscales that relate to elementary mechanisms of hydrodynamics and sediment transport(Wiggs, 2013; Courrech du Pont, 2015). Over loose granular beds, bedform emergenceis driven by a hydrodynamic instability induced by the interaction between the sand bed,flow, and sediment transport (Charru et al., 2013). On the upstream side of a bump, ero-sion takes place as the flow accelerates. Reciprocally, the flow slows down on the down-stream side where deposition occurs. However, the transition between erosion and de-position zones, associated with the location of the maximum of the sediment flux, doesnot necessarily take place at the crest of the bump. Spontaneous growth of such a bump–that is, instability–can therefore occur if its crest is located in the deposition zone (Kennedy,1963). The streamwise offset between topography and sediment flux has two contribu-tions (Andreotti et al., 2002; Kroy et al., 2002b; Fourri`ere et al., 2010; Claudin et al.,2013). First, a hydrodynamic destabilization originates from the coupling between flowinertia and dissipation, which results in a maxi- mum basal fluid shear stress located up-stream of the crest (Sykes, 1980; Hunt et al., 1988; Kroy et al., 2002a). Second, the sandflux needs a characteristic length, called the saturation length, to adapt to a spatial changein shear stress (Sauermann et al., 2001; Andreotti et al., 2010; Dur´an et al., 2011; P¨ahtzet al., 2013). This results in a stabilizing downstream lag of the maximum sand flux withrespect to the maximum of the shear stress. These balancing processes give rise to thedevelopment and propagation of sand waves at a specific wavelength and propagationspeed, associated with the most unstable mode of the pattern, with crests perpendicu-lar to the dominant wind direction.The early stage of growth and development of sedimentary ripples and dunes hasbeen theoretically studied with linear stability analyses of coupled transport and hydro-dynamic equations (Kennedy, 1963; Richards, 1980; Andreotti et al., 2002; Lagr´ee, 2003; –2–anuscript submitted to
Geophysical Research Letters
Colombini, 2004; Claudin & Andreotti, 2006; Fourri`ere et al., 2010; Devauchelle et al.,2010; Andreotti et al., 2012; Dur´an Vinent et al., 2019; Gadal et al., 2019). These anal-yses predict the incipient pattern wavelength, propagation velocity, and growth rate asfunctions of model parameters, which encode the wind and grain characteristics. For theaeolian case in particular, the dune wavelength has been shown to be proportional to thesaturation length. However, measuring the bed elevation together with sediment and fluidtransport is difficult, thereby making the direct comparison between theory and field orexperimental data rather challenging.The aerodynamic and sediment transport responses have been investigated inde-pendently of each other, and separate measurements of the saturation length and the up-wind shift of the shear stress have been carried out, in the field and in wind tunnel ex-periments (Andreotti et al., 2010; Claudin et al., 2013; Selmani et al., 2018). In contrast,few field studies addressing the early stage of aeolian dune growth are available in theliterature (Cooper, 1958; Fryberger et al., 1979; Kocurek et al., 1992; Elbelrhiti et al.,2005; Ping et al., 2014; Baddock et al., 2018). First, in situ monitoring of the evolutionof small amplitude bedforms is difficult due to the involved length and time scales (tensto hundreds of meters, days to months). Second, inherent wind variability–even in over-all unidirectional dune fields–makes application of the theory challenging. Emergenceof subaqueous sand ripples has also been experimentally investigated (Coleman & Melville,1996; Baas, 1999; Langlois & Valance, 2007; Fourri`ere et al., 2010), and more generallythe quantification of sedimentary bedforms in different environmental–including extraterrestrial–conditions in relation to hydrodynamics and sediment transport remains an active cur-rent subject of research (Lapˆotre et al., 2016; Jia et al., 2017; Lapˆotre et al., 2018; Dur´anVinent et al., 2019; Gadal et al., 2019).In this paper, we study the upwind margin of the White Sands Dune Field, wherethe dune instability leads to spatially amplifying sand waves developing downstream (Ewing& Kocurek, 2010; Phillips et al., 2019). We follow the spatiotemporal evolution of in-cipient dunes and extract their wavelength and propagation velocity, as well as their tem-poral and spatial growth rates. We then show that these four quantities all quantitativelycompare to the predictions of spatial linear stability analysis.
White Sands Dune Field is located in southern New Mexico, USA. The sand cov-ers an area of about 400 km , resulting in the largest gypsum dune field on Earth. Dom-inant winds are mainly toward the northeast and shape the sedimentary bed into trans-verse and barchan dunes, progressively turning into parabolic dunes as the vegetationcover increases (Figures 1a and 1b) (McKee, 1966; Jerolmack et al., 2012; Baitis et al.,2014). Dunes emerge on the upwind margin (Figure 1c). There, the sediment is madeof coarse, elongate, and angular grains (see Figure S2 in the supporting information),with a measured diameter d = 670 ± µ m and a bulk density ρ p = 2300 ±
100 kg m − .The saturation length, relevant in the process of dune emergence (see section 3), directlydepends on these grain properties (Andreotti et al., 2010): L sat ’ . ρ p ρ f d = 2 . ± . , (1)where ρ f ’ . − is the air mass density in ambient conditions. As one moves fur-ther into the dune field, the grain diameter and angularity both decrease, due to abra-sion and aeolian sorting (Jerolmack et al., 2011; Phillips et al., 2019). Because we re-strict our analysis of the dune development to the first kilometer along the margin, thegrain characteristics are assumed to be homogeneous. The grain roughness leads to a mea-sured avalanche slope µ = 0 . ± .
05 (see supporting information section S1).The sand flux is calculated from the hourly wind data of the weather station at Hol-loman Air Base (KHMN, 32 ◦ N, 106 ◦ W), using the method described in Courrech –3–anuscript submitted to
Geophysical Research Letters − . − . − . ◦ ]32 . . . . L a t i t ud e [ ◦ ] N Velocity at 10 m [m/s][0.0 : 2.5)[2.5 : 5.0)[5.0 : 7.5)[7.5 : 10.0)[10.0 : inf) x [km]01234 y [ k m ] N Elevation [m] . . . . . . . E l e v a t i o np r o fi l e [ m ] Lidar derived topographyMean elevation profileAnalyzed area a bcd
Figure 1.
The White Sands National Park dune field. (a) Satellite image of the White SandsDune Field (Google TM , Landsat-Copernicus). The left rose shows the wind data from 2007 to2017 (direction toward which the wind is blowing). The right rose shows the corresponding distri-bution of sand flux orientations and the resultant flux direction (brown arrow). Both agree wellwith that reported by Ewing et al. (2015); Pedersen et al. (2015). The blue area corresponds tothe location of the digital elevation measurements. (b) Digital elevation data taken in June 2007.The dashed orange line is the location of the transect shown in (c), taken along the direction ofthe resultant flux. The red star is the location of the photo shown in (d), which is a view to thesouthwest of the dune field upwind margin. –4–anuscript submitted to Geophysical Research Letters du Pont et al. (2014) (see also supporting information section S2). The wind is charac-terized by its shear velocity u ∗ , representative of the logarithmic profile inside the tur-bulent boundary layer. Its threshold value u th below which saltation cannot sustain steadytransport is estimated with u th = a p ( ρ p /ρ f ) gd = 0 .
35 m s − , and a ’ . u ∗ /u th is then about 1 . ± .
05 (all these values are gathered in Table S1). As shown in Figure 1a and documentedby Pedersen et al. (2015), the wind regime is multimodal. Southwesterly winds domi-nate the transport as noted by the nearly unimodal sand flux distribution toward thenortheast. The other modes from the north and southeast moderately impact the duneshape and migration(Swanson et al., 2016).Elevation data of the blue area in Figure 1a have been obtained using lidar-derivedtopography at five different times (June 2007, June 2008, January 2009, June 2010 andAugust 2015). Along the upwind margin, we extracted 75 dune profiles from the surfaceelevation data, spaced 50-m apart and aligned with the direction of the resulting sandflux (Figure 1b). The average topography is removed using a Butterworth high-pass fil-ter (order = 5, cutoff frequency = 0 .
005 m − ). We limit our analysis to the first dunesof the filtered profile to ensure that we stay as much as possible in the early stage of dunegrowth (red area in Figure 1c). As shown by figure 1d, these incipient dunes have verylow aspect ratios and do not exhibit any slipfaces (Phillips et al., 2019). The detrendedbed elevation exhibits a spatially amplified oscillating behaviour (Figure 2a and 2b). Wenow interpret these profiles using the theoretical framework provided by the linear sta-bility analysis developed in the next section. Here, we consider a unidirectional wind blowing at a constant shear velocity u ∗ ,over a flat sedimentary bed. Above the transport threshold velocity u th , the saturatedsediment flux induced by this flow follows a quadratic law: q sat = Ω (cid:0) u ∗ − u (cid:1) , (2)where Ω is a dimensional constant that depends on fluid and grain properties (Ungar &Haff, 1987; Iversen & Rasmussen, 1999; Creyssels et al., 2009; Dur´an et al., 2011). In nat-ural conditions, however, a sedimentary bed is never perfectly flat nor infinite, and theseirregularities can be seen as the sum of different perturbations. The purpose of the lin-ear stability analysis is precisely to study the temporal or spatial evolution of the bedin response to a perturbation of a given time or length scale. The emerging dune pat-tern is then expected to be dominated by the most unstable scale, associated with a si-nusoidal mode of wavelength λ and propagation velocity c , and whose amplitude growsin time with a rate σ or in space over a length Λ.Above such a sinusoidal bed, wind and sand flux are also modulated. As describedin the introduction, the basal shear stress is not in phase with the topography; this isquantified with two dimensionless coefficients, A and B , which represent the in-phase andin-quadrature components, respectively. They are weak (logarithmic) functions of λ (Fourri`ereet al., 2010; Charru et al., 2013) but can be, in the first approximation for our purpose,considered as constant parameters of the model. The corresponding upwind shift of thewind with respect to the bed elevation is ∼ λ B / (2 π A ). Similarly, the actual sedimentflux q is not saturated but delayed with respect to the basal shear stress, a process quan-tified by the saturation length L sat . These are the main physical mechanisms of the duneformation model from which the stability analysis is derived (see supporting informa-tion section S3 for the proper technical derivation and related theoretical figures).Consider first a large flat extent of sand. Under the action of the wind, dunes emergeeverywhere simultaneously: there is no spatial development of the pattern. A spatial si-nusoidal perturbation characterized by a given wavelength λ and an initial amplitude –5–anuscript submitted to Geophysical Research Letters − − − − B e d e l e v a t i o n [ m ] Detrended profile ± C e x/ Λ h ( x,
0) = C e x/ Λ sin(2 πx/λ ) 0 200 400 600 800Distance along profile [m] − − B e d e l e v a t i o n [ m ]
450 500 5500 . . . . . . . . . . x/λ − − h / A ∆ i x = λ ∆ i ln( H t ) = (1 / Λ) ∆ i x H it H i +1 t H it +∆ t H i +1 t +∆ t ∆ t x = c ∆ t ∆ t ln( H i ) = ς ∆ th ( x, t ) h ( x, t + ∆ t ) ± C e x/ Λ a bc Figure 2.
The spatial exponential dune growth. (a) Detrended profile corresponding to thatof Figures 1b and 1c. The black dashed lines are exponential fits to the dune crests (black dots),giving Λ = 225 m (top) and Λ = 135 m (bottom). The theoretical red profile parameters are C = 0 .
06 m, λ = 120 m, and Λ = 170 m. (b) Temporal evolution of the detrended elevationprofile, with a close-up on one crest. (c) Schematics of the spatio-temporal dune development.The theoretical profile h is defined in equation. (4).–6–anuscript submitted to Geophysical Research Letters C can grow or decay in time, at a rate σ , in response to the wind shift and the flux lag.The perturbation also propagates at a velocity c . The bed elevation along the direction x of the wind can be written as: h ( x, t ) = C e σt cos (cid:20) πλ ( x − ct ) (cid:21) . (3)Both temporal growth rate and propagation speed can be computed as functions of λ from the analysis of the equations coupling the flow, sediment transport and bed evo-lution, constituting the dispersion relation of sand waves (supporting information sec-tion S3). Positive values of the growth rate are associated with unstable perturbations,and these are typically for large values of λ . Conversely, small wavelengths are stable,with σ < x = 0. Dunes emerge from the selective am-plification of disturbances propagating downwind from the field entrance. This resultsin a spatial development of the pattern. There is no temporal growth: at a given loca-tion, the amplitude of the bed oscillation is the same at any time. The form of a sinu-soidal mode of initial amplitude C can be written as: h ( x, t ) = C e x/ Λ cos (cid:20) πλ ( x − ct ) (cid:21) , (4)where Λ − is the spatial growth rate of the dunes.Neutral modes are the same in both spatial and temporal analysis. Their wavelength λ c is characterised by vanishing growth rates σ ( λ c ) = 0 and Λ − ( λ c ) = 0, such that, λ c = 2 π AB − µ (cid:18) u th u ∗ (cid:19) L sat , (5)and separates growing ( λ > λ c ) from decaying ( λ < λ c ) perturbations. It can thusbe interpreted as a minimal dune size.Performing the temporal linear stability analysis (denoted by subscript T), in thelimit L sat /λ c (cid:28)
1, the characteristics of the fastest growing perturbation read: λ T ∼ λ c , (6) σ T ∼ QL A (cid:18) πL sat λ T (cid:19) , (7) c T ∼ QL sat A πL sat λ T , (8)where Q = Ω u ∗ gives the characteristic scale of the sediment flux associated with thewind regime (Fourri`ere et al., 2010; Gadal et al., 2019).Conversely, spatial growth rate reaches a maximum 1 / Λ S at some specific valueof the wavelength, noted λ S , corresponding to the propagation velocity c S . Unfortunately,no simple analytical and accurate formulae like (6-8) can be derived for these quantities(see supporting information section S3). Temporal and spatial analyses are consistent,and we typically find λ S ’ . λ T and c S ’ . c T . The numerical factors in these re-lations do not vary by more than a few percent upon changing the parameters A and B within a reasonable range of values. In the spatial development of the instability, anindividual bump grows in height while propagating downwind at a constant velocity. Itsamplitude therefore varies exponentially with respect to time, and one can define a tem-poral growth rate ς S ≡ c S / Λ S (inset of Figure 2c) in the frame of reference of the bumpi.e. a Lagrangian growth rate. The theoretical analysis provides the approximate rela-tion ς S ’ . σ T . Both ς S and Λ S can be measured from the field data, providing equiv-alent information. –7–anuscript submitted to Geophysical Research Letters
After removing the average topography (Figure 1c), the detrended bed elevationprofile exhibits an exponentially amplifying sinusoidal shape, as predicted by the spa-tial linear theory for dune emergence (Figure 2a). Using these profiles and their tempo-ral evolution (Figure 2b), we have extracted the three independent characteristics of thepattern ( λ , c and Λ or ς ) using two different methods. For each profile, we either lookat each peak separately (peak-to-peak method), or extract quantities averaged over thewhole profile (global approach).The wavelength λ is computed by autocorrelation of the bed elevation profile (globalmethod), and from the spacing between two adjacent peaks (peak-to-peak method). Thefit of an exponential to the peaks of each profile gives the spatial growth length Λ (globalmethod, see dashed line in Figure 2a). The spatial growth length is alternatively com-puted from the difference in height between two adjacent peaks (peak-to-peak method,see Figure 2c).The Lagrangian growth rate ς and propagation velocity c are obtained by fittingexponential and linear functions to the temporal variation of the dune height and po-sition, respectively (see equations 43 and 44 of the supporting information section S3).The peak-to-peak method looks at the height and position of each peak separately (seeinset of Figure 2c). For a global measurement, the average propagation speed can be de-termined from the cross-correlation curve between the same profiles at different times,and the average Lagrangian growth rate from the temporal evolution of the bed eleva-tion standard deviation.Importantly, to extract the average values of ς and c , we take into account the tem-poral variations of the characteristic sand flux Q . Indeed, both quantities vary in timeproportionally to Q . We also remove periods of time when the wind is below the trans-port threshold (see supporting information section S4). The output of the analysis of the 75 transects is shown in Figure 3. Both peak-to-peak and global methods exhibit similar distributions for the wavelength λ , the prop-agation velocity c , the Lagrangian growth rate ς and the growth length Λ, with clear dom-inant (most probable) values. The incipient dune wavelength and growth length are bothon the order of a hundred meters; their propagation velocity is around 5 m yr − and theirgrowth rate is about 0 .
015 yr − . These values, as well as the typical dispersion aroundthem (i.e. the width of these distributions), are more precisely reported in Table S2. For λ and c , our results are consistent with the measurements of Phillips et al. (2019), madeon a single elevation profile. Their dispersion is on the order of 20%, because these quan-tities can be measured with a good accuracy, especially with the global method usingcorrelation. As the measurement of Λ and ς is more delicate, the corresponding distri-butions are more dispersed. The peak-to-peak method is actually sensitive to the behaviourof individual peaks, that can respond to various types of local disturbances. For exam-ple, they may induce irregularities in the spacing between the peaks, or asymmetry be-tween positive and negative detrended topography (Figure 2b). As a result, a few neg-ative values of the velocity, growth rate and characteristic growth length are reported.Nevertheless, these data provide reliable and meaningful statistics to test the theory, whichmust be able to account for those four quantities concomitantly. The free parameters ofthe linear analysis are the hydrodynamic coefficients A and B , as the others are set in-dependently with the wind and sediment properties (see Table S1).The incipient dune wavelength peaks around 120 m, and is therefore significantlylarger than the usually reported value ( ’
20 m, for sand particles of size 180 µ m) (Elbelrhitiet al., 2005). How can the theory reproduce such a large value? First, as the wavelength –8–anuscript submitted to Geophysical Research Letters is proportional to the saturation length, and thus to the grain size (equations 1, 5 and 6),the presence here of much coarser grains provides a factor 670 / ’ . u ∗ /u th ’ .
3, we can expect an additional factor of 2 with respect to a situationfar above threshold. The simultaneous fit of the four quantities predicted by the spatiallinear stability analysis to the data allows us to reproduce quantitatively the dominantwavelength and growth rate. The predicted values fall in the peak of the distributions(red lines in Figures 3a and 3c). This adjustment however overestimates the growth lengthand the propagation velocity, whose predictions exceed the dominant values by an amountcomparable to the peaks’ width (red lines in Figures 3b and 3d).This discrepancy can be understood by questioning our approximation of a uni-directional wind. A finer analysis of the flux rose shows in fact secondary winds, withnon-zero components perpendicular to the crest toward the southwest. Reversing windshave cumulative effects on the growth rate and the selected wavelength. They howeverpartially cancel each other out, impeding the propagation, and thus the spatial devel-opment of the dune pattern. Such a process is supported by observations of reactivationsurfaces formed by reversing winds in the stratigraphy at White Sands (Phillips et al.,2019; Kocurek et al., 2007). The value of the characteristic sand flux Q (given in Table S1)is a time average that does not account for changes in wind orientation. The ratio be-tween the scalar and vector averages of Q , taking into account the variations in orien-tation of the sand fluxes over time, is about 0 . The exponentially amplified sinusoidal behaviour of the White Sands bed eleva-tion profiles is field evidence for the spatial development of the dune instability. The lin-ear analysis is able to quantitatively reproduce the three characteristics of the emergentpattern: dune wavelength, propagation velocity and growth length (or equivalently, growthrate). To obtain this agreement, the two hydrodynamic coefficients where adjusted, re-sulting in A = 3 . ± . B = 1 . ± .
3. All the other parameters of the theorywere fixed independently from sediment and wind data, either by direct measurement(grain diameter and density, avalanche slope) or using well calibrated relationships (sat-uration length, sediment flux). Note that the uncertainty on the determined coefficientsis dominated by the those of the fixed parameters, rather than by the dispersion in themeasurements of c , λ and Λ. Overall, the value of L sat , involved in λ and Λ, directly af-fects the estimate of B / A , while the value of Q , involved in c , mostly affects A .Importantly, the concomitant agreement of c , λ and Λ is a stringent test of the the-ory. This study is therefore a step forward in the general ‘dune inverse problem’, try-ing to infer, for example, grain or wind properties from dune characteristics (Fenton etal., 2014a, 2014b; Ewing et al., 2015; Runyon et al., 2017; Fernandez-Cascales et al., 2018).It is remarkable that the resulting values determined from our field data are very closeto those directly measured by Claudin et al. (2013) on a single 40 m long dune ( A ’ . B ’ . –9–anuscript submitted to Geophysical Research Letters
50 100 150 200Wavelength λ [m]0 . . . . . . . P D F × − −
200 0 200 400 600Growth length Λ [m]02468 P D F × − − . − . . . ς [yr − ]051015 P D F c [m yr − ]01234 P D F × − a bc d Figure 3.
Distributions of incipient dune time and length scales. Blue and green distributionsshows the results of the global and peak-to-peak methods, respectively. Errorbars gives the rangeof values obtained from the spatial linear stability analysis with A = 3 . ± . B = 1 . ± . λ c is in grey. –10–anuscript submitted to Geophysical Research Letters
A limitation of the linear theory is of course the presence of non-linear effects. Theyoccur when the aspect ratio of the sand waves becomes too large or when the dunes in-teract with each other, so that each bed perturbation cannot be considered as indepen-dent of the others. Bumps with aspect ratios above ’ /
13 are expected to start to de-velop flow recirculation on their downwind side, usually associated with the formationof an avalanche slip face (Fourri`ere et al., 2010). In the region we have studied, the wavesfurthest inside the dune field could reach aspect ratios of about 1 /
10, but no slip faceswere observed. Similar to Phillips et al. (2019), we also recognize the coexistence of mul-tiple wavelengths at the upwind side of the profiles (associated with different celeritiesand growth rates or lengths), and these are partly the cause of the distribution widthsin Figure 3. We could not, however, infer from these data signs of interactions, such ascollisions, coalescence or ejection (Hersen & Douady, 2005; Katsuki et al., 2005; Gao etal., 2015b; Bacik et al., 2020).Although studied here at the boundary of a dune field, the spatial development ofthe dune instability is also present on pre-existing large dunes, providing a similar up-wind boundary condition in terms of sand availability. As a matter of fact, the emergenceof bed oscillations on the flanks of barchans has been proposed as a key mechanism tounderstand their stability, as these superimposed waves eventually grow until they canbreak from the horns, causing large sand losses (Elbelrhiti et al., 2005; Zhang et al., 2010;Lee et al., 2019). Likewise, in narrow bidirectional wind regimes, the growth of the in-stability over elongating linear dunes breaks them into trains of barchans (Gao et al.,2015). This work therefore provides a reliable base to study the stability of large dunesand thus the formation of large-scale structures inside dune fields (Worman et al., 2013;Gadal et al., 2020).
Acknowledgments
This collaborative work was initiated at the International Conference on Aeolian Research(ICAR), 2018, Bordeaux, France. We acknowledge financial support from the Univ-EarthSLabEx program of Sorbonne Paris Cit´e (grants ANR-10-LABX-0023 and ANR-11-IDEX-0005-02) and the French National Research Agency (grant ANR-17-CE01-0014/SONO).Cl´ement Narteau acknowledges support from the National Science Center of Poland (grant2016/23/B/ST10/01700). Philippe Claudin acknowledges support from TOAD (The Orig-ine of Aeolian Dunes) project as an external partner (Natural Environment Research Coun-cil, UK and National Science Foundation, USA; NE/R010196NSFGEO-NERC). RyanC. Ewing acknowledges support from the White Sands National Monument through NPS-GC-CESU Cooperative Agreement https://mesowest.utah.edu/cgi-bin/droman/mesomap.cgi?state=NM&rawsflag=3 ). Topographic data can be foundin the public repositories Texas Data Repository and OpenTopography (2007-doi:10.18738/T8/WUNF0G,2008-doi:10.18738/T8/HQVUSX, 2009-doi:10.5069/G9Q23X5P, 2010-doi:10.5069/G97D2S2D,2015: https://portal.opentopography.org/usgsDataset?dsid=USGS LPC NM WhiteSands2015 LAS 2017 ). Supplementary figures, tables, and text can be found in the support-ing information. We thank the anonymous reviewers for their careful reading of our manuscriptand their insightful comments and suggestions.
References
Andreotti, B., Claudin, P., Devauchelle, O., Dur´an, O., & Fourri`ere, A. (2012). Bed-forms in a turbulent stream: ripples, chevrons and antidunes.
Journal of FluidMechanics , , 94–128.Andreotti, B., Claudin, P., & Douady, S. (2002). Selection of dune shapes and veloc-ities: Part 2: A two-dimensional modelling. European Physical Journal B , ,341–352.Andreotti, B., Claudin, P., & Pouliquen, O. (2010). Measurements of the aeolian –11–anuscript submitted to Geophysical Research Letters sand transport saturation length.
Geomorphology , , 343–348.Baas, J. H. (1999). An empirical model for the development and the equilibriummorphology of current ripples in fine sand. Sedimentology , , 123–138.Bacik, K. A., Lovett, S., Colm-cille, P. C., & Vriend, N. M. (2020). Wake inducedlong range repulsion of aqueous dunes. Physical Review Letters , (5),054501.Baddock, M. C., Nield, J. M., & Wiggs, G. F. S. (2018). Early-stage aeolian pro-todunes: Bedform development and sand transport dynamics. Earth SurfaceProcesses and Landforms , , 339–346.Baitis, E., Kocurek, G., Smith, V., Mohrig, D., Ewing, R. C., & Peyret, A.-P. B.(2014). Definition and origin of the dune-field pattern at White Sands, NewMexico. Aeolian Research , , 269–287.Charru, F., Andreotti, B., & Claudin, P. (2013). Sand ripples and dunes. AnnualReview of Fluid Mechanics , , 469–493.Claudin, P., & Andreotti, B. (2006). A scaling law for aeolian dunes on Mars, Venus,Earth, and for subaqueous ripples. Earth and Planetary Science Letters , ,30–44.Claudin, P., Wiggs, G. F. S., & Andreotti, B. (2013). Field evidence for the upwindvelocity shift at the crest of low dunes. Boundary-layer Meteorology , , 195–206.Coleman, S. E., & Melville, B. W. (1996). Initiation of bed forms on a flat sand bed. Journal Hydraulic Engineering , , 301–310.Colombini, M. (2004). Revisiting the linear theory of sand dune formation. Journalof Fluid Mechanics , , 1–16.Cooper, W. S. (1958). Coastal sand dunes of oregon and washington (Vol. 72). Geo-logical Society of America.Courrech du Pont, S. (2015). Dune morphodynamics.
Compte-Rendus de Physique , , 118–138.Courrech du Pont, S., Narteau, C., & Gao, X. (2014). Two modes for dune orienta-tion. Geology , , 743–746.Creyssels, M., Dupont, P., Ould el Moctar, A., Valance, A., Cantat, I., Jenkins,J. T., . . . Rasmussen, K. R. (2009). Saltating particles in a turbulent bound-ary layer: experiment and theory. Journal of Fluid Mechanics , , 47–74.Devauchelle, O., Malverti, L., Lajeunesse, E., Lagr´ee, P.-Y., Josserand, C., &Nguyen Thu-Lam, K.-D. (2010). Stability of bedforms in laminar flows withfree surface: from bars to ripples. Journal of Fluid Mechanics , , 329–348.Dur´an, O., Claudin, P., & Andreotti, B. (2011). On aeolian transport: Grain-scaleinteractions, dynamical mechanisms and scaling laws. Aeolian Research , ,243–270.Dur´an Vinent, O., Andreotti, B., Claudin, P., & Winter, C. (2019). A unified modelof ripples and dunes in water and planetary environments. Nature Geoscience , , 345–350.Elbelrhiti, H., Claudin, P., & Andreotti, B. (2005). Field evidence for surface-wave-induced instability of sand dunes. Nature , , 720–723.Ewing, R. C., & Kocurek, G. (2010). Aeolian dune interactions and dune?fieldpattern formation: White Sands Dune Field, New Mexico. Sedimentology , ,1199–1219.Ewing, R. C., McDonald, G. D., & Hayes, A. G. (2015). Multi-spatial analysis of ae-olian dune-field patterns. Geomorphology , , 44–53.Fenton, L. K., Michaels, T. I., & Beyer, R. A. (2014a). Inverse maximum grossbedform-normal transport 1: How to determine a dune-constructing windregime using only imagery. Icarus , , 5–14.Fenton, L. K., Michaels, T. I., Chojnacki, M., & Beyer, R. A. (2014b). Inverse maxi-mum gross bedform-normal transport 2: Application to a dune field in GangesChasma, Mars and comparison with HiRISE repeat imagery and MRAMS. –12–anuscript submitted to Geophysical Research Letters
Icarus , , 47–63.Fernandez-Cascales, L., Lucas, A., Rodriguez, S., Gao, X., Spiga, A., & Narteau, C.(2018). First quantification of relationship between dune orientation and sedi-ment availability, olympia undae, mars. Earth and Planetary Science Letters , , 241–250.Fourri`ere, A., Claudin, P., & Andreotti, B. (2010). Bedforms in a turbulent stream:formation of ripples by primary linear instability and of dunes by nonlinearpattern coarsening. Journal of Fluid Mechanics , , 287–328.Fryberger, S. G., Ahlbrandt, T. S., & Andrews, S. (1979). Origin, sedimentaryfeatures, and significance of low-angle eolian ”sand sheet” deposits, Great SandDunes National Monument and vicinity, Colorado. Journal of SedimentaryResearch , (3), 733–746.Gadal, C., Narteau, C., Courrech du Pont, S., Rozier, O., & Claudin, P. (2019). In-cipient bedforms in a bidirectional wind regime. Journal of Fluid Mechanics , , 490–516.Gadal, C., Narteau, C., Courrech du Pont, S., Rozier, O., & Claudin, P. (2020). Pe-riodicity in fields of elongating dunes. Geology , (4), 343–347.Gao, X., Narteau, C., & Rozier, O. (2015b). Development and steady states oftransverse dunes: A numerical analysis of dune pattern coarsening and giantdunes. Journal of Geophysical Research: Earth Surface , (10), 2200–2219.Gao, X., Narteau, C., Rozier, O., & Courrech Du Pont, S. (2015). Phase diagrams ofdune shape and orientation depending on sand availability. Scientific reports , , 1–12.Hersen, P., & Douady, S. (2005). Collision of barchan dunes as a mechanism of sizeregulation. Geophysical Research Letters , (21).Hunt, J. C. R., Leibovich, S., & Richards, K. J. (1988). Turbulent shear flowsover low hills. Quarterly Journal of the Royal Meteorological Society , ,1435–1470.Iversen, J. D., & Rasmussen, K. R. (1999). The effect of wind speed and bed slopeon sand transport. Sedimentology , , 723–731.Jerolmack, D. J., Ewing, R. C., Falcini, F., Martin, R. L., Masteller, C., Phillips, C.,. . . Buynevich, I. (2012). Internal boundary layer model for the evolution ofdesert dune fields. Nature Geoscience , (3), 206–209.Jerolmack, D. J., Reitz, M. D., & Martin, R. L. (2011). Sorting out abrasion in agypsum dune field. Journal of Geophysical Research , , F02003.Jia, P., Andreotti, B., & Claudin, P. (2017). Giant ripples on comet 67p/churyumov-gerasimenko sculpted by sunset thermal wind. Proceedings of the NationalAcademy of Sciences USA , , 2509–2514.Katsuki, A., Nishimori, H., Endo, N., & Taniguchi, K. (2005). Collision dynamics oftwo barchan dunes simulated using a simple model. Journal of the Physical So-ciety of Japan , (2), 538–541.Kennedy, J. F. (1963). The mechanics of dunes and antidunes in erodible-bed chan-nels. Journal of Fluid Mechanics , , 521–544.Kocurek, G., Carr, M., Ewing, R. C., Havholm, K. G., Nagar, Y. C., & Singhvi,A. K. (2007). White Sands dune field, New Mexico: age, dune dynamics andrecent accumulations. Sedimentary Geology , , 313–331.Kocurek, G., Townsley, M., Yeh, E., Havholm, K. G., & Sweet, M. L. (1992). Duneand dune-field development on Padre Island, Texas, with implications forinterdune deposition and water-table-controlled accumulation. Journal ofSedimentary Research , (4), 622–635.Kroy, K., Sauermann, G., & Herrmann, H. J. (2002a). Minimal model for aeoliansand dunes. Physical Review E , , 031302.Kroy, K., Sauermann, G., & Herrmann, H. J. (2002b). Minimal Model for SandDunes. Physical Review Letters , , 054301. –13–anuscript submitted to Geophysical Research Letters
Lagr´ee, P.-Y. (2003). A triple deck model of ripple formation and evolution.
Physicsof Fluids , , 2355–2368.Langlois, V., & Valance, A. (2007). Formation and evolution of current ripples ona flat sand bed under turbulent water flow. European Physical Journal E , ,201–208.Lapˆotre, M. G. A., Ewing, R. C., Lamb, M. P., Fischer, W. W., Grotzinger, J. P.,Rubin, D. M., . . . others (2016). Large wind ripples on Mars: A record ofatmospheric evolution. Science , (6294), 55–58.Lapˆotre, M. G. A., Ewing, R. C., Weitz, C. M., Lewis, K. W., Lamb, M. P.,Ehlmann, B. L., & Rubin, D. M. (2018). Morphologic diversity of Martianripples: implications for large-ripple formation. Geophysical Research Letters , , 10229–10239.Lee, D. B., Ferdowsi, B., & Jerolmack, D. J. (2019). The imprint of vegetation ondesert dune dynamics. Geophysical Research Letters , (21), 12041–12048.McKee, E. D. (1966). Structures of dunes at White Sands National Monument, NewMexico (and a comparison with structures of dunes from other selected areas)1. Sedimentology , (1), 3–69.Narteau, C., Zhang, D., Rozier, O., & Claudin, P. (2009). Setting the length andtime scales of a cellular automaton dune model from the analysis of superim-posed bed forms. Journal of Geophysical Research , , F03006.P¨ahtz, T., Kok, J. F., Parteli, E. J. R., & Herrmann, H. J. (2013). Flux saturationlength of sediment transport. Physical Review Letters , , 218002.Pedersen, A., Kocurek, G., Mohrig, D., & Smith, V. (2015). Dune deformation ina multi-directional wind regime: White Sands Dune Field, New Mexico. EarthSurface Processes and Landforms , , 925–941.Phillips, J. D., Ewing, R. C., Bowling, R., Weymer, B. A., Barrineau, P., Nittrouer,J. A., & Everett, M. E. (2019). Low-angle eolian deposits formed by protodunemigration, and insights into slipface development at White Sands dune field,New Mexico. Aeolian Research , , 9–26.Ping, L., Narteau, C., Dong, Z., Zhang, Z., & Courrech du Pont, S. (2014). Emer-gence of oblique dunes in a landscape-scale experiment. Nature Geoscience , ,99–103.Richards, K. J. (1980). The formation of ripples and dunes on an erodible bed. Journal of Fluid Mechanics , , 597–618.Runyon, K., Bridges, N., Ayoub, F., Newman, C., & Quade, J. (2017). An inte-grated model for dune morphology and sand fluxes on Mars. Earth and Plane-tary Science Letters , , 204–212.Sauermann, G., Kroy, K., & Herrmann, H. J. (2001). Continuum saltation model forsand dunes. Physical Review E , , 031305.Selmani, H., Valance, A., Ould el Moctar, A., Dupont, P., & Zegadi, R. (2018). Ae-olian sand transport in out-of-equilibrium regimes. Geophysical Research Let-ters , , 1838–1844.Swanson, T., Mohrig, D., & Kocurek, G. (2016). Aeolian dune sediment flux vari-ability over an annual cycle of wind. Sedimentology , , 1753–1764.Sykes, R. I. (1980). An asymptotic theory of incompressible turbulent boundarylayer flow over a small hump. Journal of Fluid Mechanics , , 647–670.Ungar, J. E., & Haff, P. K. (1987). Steady state saltation in air. Sedimentology , ,289–299.Wiggs, G. F. S. (2013). Dune morphology and dynamics. Treatise on Geomorphol-ogy , , 201–218.Worman, S. L., Murray, A. B., Littlewood, R., Andreotti, B., & Claudin, P. (2013).Modeling emergent large-scale structures of barchan dune fields. Geology , ,1059–1062.Zhang, D., Narteau, C., & Rozier, O. (2010). Morphodynamics of barchan andtransverse dunes using a cellular automaton model. Journal of Geophysical Re- –14–anuscript submitted to
Geophysical Research Letters search: Earth Surface , (F3). –15– EOPHYSICAL RESEARCH LETTERS
Supporting Information for “Spatial and temporaldevelopment of the dune instability”
C. Gadal , C. Narteau , R.C. Ewing , A. Gunn , D. Jerolmack , B.Andreotti , P. Claudin Institut de physique du globe de Paris, Universit´e de Paris, CNRS, Paris, France Department of Geology and Geophysics, Texas A&M University, College Station, USA Department of Earth and Environmental Science, University of Pennsylvania, Philadelphia, USA Laboratoire de Physique, ENS - PSL Research University, Universit´e de Paris, CNRS, Sorbonne Universit´e, Paris, France Physique et M´ecanique des Milieux H´et´erog`enes, ESPCI Paris - PSL Research University, Universit´e de Paris, CNRS, SorbonneUniversit´e, Paris, France
Corresponding author: C. Gadal, Institut de Physique du globe de Paris, 1 rue Jussieu, Paris,France ([email protected])
October 5, 2020, 2:44am a r X i v : . [ phy s i c s . g e o - ph ] O c t - 2 GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
Contents of this file • Introduction • Data availability statement • Text S1 to S5:1. Text S1. Measurement of the avalanche angle2. Text S2. Wind data processing3. Text S3. Linear stability analysis4. Text S4. Temporal fit of the measured quantities5. Text S5. Distributions of the hydrodynamic coefficients and uncertainties • Figures S1 to S6 • Tables S1 to S2
October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 3
Introduction
In this supporting information, we provide details on data measurement and processingfor the avalanche angle (Text S1) and wind data (Text S2). We also derive the temporaland spatial linear stability analysis (Text S4). The spatial analysis technically differs fromthe standard temporal analysis and has never been published before in the literature.We also explain how to appropriately scale time when comparing data and theory whensediment transport is intermittent (Text S4), and how the hydrodynamic parameters inthe analysis can be estimated from the data (Text S5). This document is prefaced withinformation on where to find the data, and ends with supporting figures pertaining toText S1-S5 (Figures S1-S6) and tables summarizing all numerical quantities useful for ouranalysis (Tables S1-S2).
Data availability statement
The wind data is analyzed for the period 2007-2010, based on the data collectedby the United States Air Force, and hosted by MesoWest under the station code‘KHNM’ here: https://mesowest.utah.edu/cgi-bin/droman/mesomap.cgi?state=NM&rawsflag=3 . All of the avalanche slope data is reported in Figure S1. The topo-graphic data information is as follows: • June 2007 (DOI:10.18738/T8/WUNF0G) hosted by Texas Data Repository, collectedby Gary Kocurek, and funded by the National Park Service. • June 2008 (DOI:10.18738/T8/HQVUSX) hosted by Texas Data Repository, collectedby Gary Kocurek, and funded by the National Park Service.
October 5, 2020, 2:44am - 4
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY • January 2009 (DOI:10.5069/G9Q23X5P) hosted by OpenTopography, collected byGary Kocurek, and funded by the National Center for Airborne Laser Mapping. • June 2010 (DOI:10.5069/G97D2S2D) hosted by OpenTopography, collected by RyanEwing, and funded by the National Center for Airborne Laser Mapping. • August 2015 is hosted by OpenTopography and collected as part of the UnitedStates Geological Survey 3D Elevation Program. Available here: https://portal.opentopography.org/usgsDataset?dsid=USGS LPC NM WhiteSands 2015 LAS 2017 . Text S1. Measurement of the avalanche slope µ The avalanche slope µ of a granular material is measured from the shape of a conicalpile built with these grains. Its value slightly depends on the grains’ shape, and, because itenters the parameters of the linear stability analysis, we have measured it for the gypsumparticles that compose the dunes at the the upwind margin of the White Sands. Theexperimental set-up is presented in Figure S1. We fill a cylindrical tube with the grainsand slowly pull it up. The grains flow out of the tube bottom, forming, in a quasi-staticway, a conical pile. The slope of the pile evolves through a succession of avalanches. Wetake pictures after every incremental lift of the tube. We have so obtained data for 39cones of different volumes. On these pictures, we detect the edges of the pile by simpleimage processing, and fit them by straight lines to measure the corresponding slope. Theresults are shown in Figure S1. The values are nicely independent of the pile height, as itshould for homothetic cones, showing that, finite size effects (the piles are small) can beneglected. These data allow us to obtain a mean value: µ = 0 . ± . , (1) October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 5 corresponding to an avalanche angle ’ ◦ . It is slightly larger than the usual value of35 ◦ , in agreement with the elongated and angular shape of the grains (see figure S2). Thedispersion of the data relates to the metastability of the avalanching process. Text S2. Wind data processing
The sand fluxes are derived from the wind data of the KHMN weather station at theHolloman air base (32 ◦ N, 106 ◦ W) corresponding to the topography measurements,between 2007 and 2010. The method of analysis used here is based on that described inthe supplementary material of Courrech du Pont, Narteau, and Gao (2014). The dataare hourly measurements of the average wind velocity u t and direction θ t at 10 m height,with velocity and angle bins of 0 . ◦ , respectively. The sub-letter t indicateshere the time of the corresponding data.The sediment flux depends on the shear velocity u ∗ , characteristic to the logarithmicvelocity vertical profile of the wind u ( z ) = u ∗ /κ ln ( z/z ) inside the turbulent boundarylayer. κ = 0 . z is the aerodynamic roughness. At time t , this law of the wall is reciprocally used to estimate the shear velocity from the windmeasurement as: u ∗ ,t = κ u t ln (cid:16) zz (cid:17) , (2)where z is set to the height at which the wind data have been recorded. Here we take z ’ u th to maintain steady salta-tion. Its value is empirically given by u th = 0 . r ρ p ρ f gd, (3) October 5, 2020, 2:44am - 6
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY where ρ p and ρ f are the grain and air density, d the grain diameter and g the gravitationalacceleration (Iversen & Rasmussen, 1999). The saturated flux can then be computed ateach time step using the transport law (Ungar & Haff, 1987): q sat ,t = ( Ω (cid:0) u ∗ ,t − u (cid:1) e θ t if u ∗ ,t > u th , , (4)where e θ t = (cos θ t , sin θ t ) is a unit vector, and Ω = 25 ( ρ f /ρ p ) p d/g is a dimensionalconstant calibrated by measurements and numerical simulations (see review of Dur´an,Claudin, and Andreotti (2011)). From this time series, one can define the resultant sandflux. Likewise, we also compute the characteristic sand flux Q as: Q t = ( Ω u ∗ ,t e θ t if u ∗ ,t > u th , . (5)Finally, as the sand flux direction is mostly unimodal, dune growth and migration mainlyresult from the component of the fluxes aligned with the resultant sand flux direction β ,perpendicular to the dune crests. The relevant characteristic sand flux for the study ofthe spatial development of the dune instability then become: Q t = ( Ω u ∗ ,t cos ( θ t − β ) if u ∗ ,t > u th , . (6)Note that Q t can have positive or negative values depending on the orientation of the sandflux at time with respect to the resultant (average) sand flux direction. Quantification ofthe impact of the variations in wind direction can then be done by looking at the ratiobetween scalar and vector average of the characteristic sand flux. Using the wind data weobtain h| Q t |i = 38 ± yr − , (7) h Q t i = 23 ± yr − , (8) October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 7 resulting in h Q t ih| Q t |i ’ . . (9) Text S3. Linear stability analysis
We derive here the linear stability analysis of a flat sedimentary bed sheared by aconstant fluid flow, detailing the differences between temporal and spatial analyses. Thecalculations closely follow the framework presented in the review of Charru, Andreotti,and Claudin (2013), leading to the characteristics of the most unstable mode.As all unstable modes are found to be propagating downwind, the instability is saidto be convective : it grows while propagating. Then, only for an artificial situation withperiodic boundary conditions, for example in numerical simulations, the instability cansimultaneous grow in time, everywhere. Likewise, if one could flatten the sand bed at akilometer scale, one would see the spatially-homogeneous emergence of the instability intime. Nevertheless, the dunes generated from noise at the upstream boundary conditionwould eventually invade, by propagation, a larger and larger domain.
General framework and notations — The evolution of the bed elevation profile h ( x, t )is governed by the mass conservation (Exner) equation ∂ t h + ∂ x q = 0 , (10)where q is the sediment flux. In the steady homogenous case, the sediment flux takes itssaturated value, which is a linear function of the basal shear stress τ b ≡ ρ f u ∗ q sat = Ω ρ f ( τ b − τ th ) , (11) October 5, 2020, 2:44am - 8
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY where τ th ≡ ρ f u is the threshold value of the stress, below which transport vanishes.In non-homogeneous situations, the flux lags behind its saturated value, following therelaxation equation L sat ∂ x q = q sat − q, (12)where L sat is the saturation length.In the context of the linear stability analysis, we consider a bed perturbation of theform h ( x, t ) = He i ( kx − ωt ) , (13)where where k is the wavenumber and ω the associated pulsation. H is the amplitudeof the perturbation, and we always assume to work in the limit kH (cid:28)
1. On such asinusoidal bed, the associated basal shear perturbation writes in Fourier notation:ˆ τ b = ρ f u ∗ ( A + i B ) kH. (14) A and B are hydrodynamic coefficients representing the in-phase and in-quadrature partof the shear stress with respect to the topography. Grouped together, the above equationslead to the following dispersion relation: L Q ω = ( L sat k ) iL sat k ( A + i B ) , (15)where Q = Ω u ∗ is the characteristic sand flux. If one accounts for the dependence of thethreshold on the bed slope ∂ x h , one should simply replace the coefficient B by: B µ = B − µ (cid:18) u th u ∗ (cid:19) , (16)where µ is the avalanche slope. October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 9
Generally speaking, the wavenumber and pulsation are complex numbers: k = k r + ik i , (17) ω = ω r + iω i , (18)whose real parts ( k r , ω r ) represent spatial and temporal oscillations, and whose imaginaryparts ( k i , ω i ) represent growth or decay of the amplitude of the perturbation. For sakeof clarity, the pulsation and wavenumber will thereby only refer to their real parts, whileone defines the temporal and spatial growth rates from their imaginary parts. From nowon, the notations will then be the following: k r ≡ k, (19) ω r ≡ ω, (20) ω i ≡ σ, (21) − k i ≡ / Λ , (22)where Λ is the characteristic growth length. We also define the propagation velocity ofthe pattern as c = ω/k . Temporal linear stability analysis — We first analyse the temporal stability of system.We consider a spatially infinite system, where we disturb the bed at a time t = 0 by aspatial oscillation of given wavelength λ = 2 π/k . The system remains spatially homoge-neous, which implies 1 / Λ = 0: there is neither spatial growth nor decay. We thereforelook for the temporal response of the system to this perturbation, to compute its growth
October 5, 2020, 2:44am - 10
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY rate. The real and imaginary parts of the dispersion relation (15) can be expressed as: L sat Q c ( k ) = L sat k L sat k ) ( A + B µ L sat k ) , (23) L Q σ ( k ) = ( L sat k ) L sat k ) ( B µ − A L sat k ) . (24)Now we find the most unstable wavenumber k T (Figure S3(a)) corresponding to the wave-length λ T = 2 π/k T at the maximum σ , assuming that A and B are constant (they arein fact weak (logarithmic) functions of k in the rough turbulent regime of interest here(Fourri`ere et al., 2010)). This most unstable wavelength writes in the limit L sat k (cid:28) λ T ∼ π AB µ L sat . (25)The corresponding growth rate and propagation velocity are: σ T = σ ( k T ) ∼ QL A (cid:18) πλ T L sat (cid:19) = QL B µ A , (26) c T = c ( k T ) ∼ QL sat A (cid:18) πλ T L sat (cid:19) = QL sat B µ . (27)The Lagrangian growth length can be defined as c T /σ T . Spatial linear stability analysis — Likewise, one can investigate the spatial stability ofsystem. We perturb in time a specific point x = 0 of the sediment bed at a given pulsation ω . As the system is temporally homogeneous at the beginning, so is its response, whichimplies σ = 0 (no temporal growth or decay). We therefore look for the spatial responseof the system to this perturbation, to compute its growth length. The dispersion relation(15) can be rewritten as:( A + i B µ ) L (cid:0) k − i Λ − (cid:1) − i L Q ωL sat (cid:0) k − i Λ − (cid:1) − L Q ω = 0 , (28) October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 11 which solves into: L sat (cid:0) k − i Λ − (cid:1) = 12 ( A + i B µ ) " i L Q ω ± s L Q ω (cid:20) A + i B µ ) − L Q ω (cid:21) . (29)The two roots k + and k − relate to waves propagating downstream and upstream, respec-tively. Splitting the above expression into real and imaginary part, we obtain: L sat k ( ω ) = 12 (cid:0) A + B µ (cid:1) (cid:18) B µ L Q ω ± F ( ω ) [ AG ( ω ) + B µ H ( ω )] (cid:19) , (30) L sat Λ − ( ω ) = 12 (cid:0) A + B µ (cid:1) (cid:18) −A L Q ω ± F ( ω ) [ B µ G ( ω ) − AH ( ω )] (cid:19) , (31)where F ( ω ) = s L Q ω (cid:18) A − L Q ω (cid:19) , (32) G ( ω ) = vuuuuuut vuuuuut B µ A − L Q ω , (33) H ( ω ) = vuuuuuutvuuuuut B µ A − L Q ω − . (34)The spatial growth rate 1 / Λ is plotted for both branches in Figure S3b as a function of thepulsation ω . Only the positive branch k + corresponding to waves propagating downstreamexhibits positive values, and can grow. It also shows a single maxima, corresponding tothe most unstable pulsation ω S . The corresponding wavelength, propagation velocity and October 5, 2020, 2:44am - 12
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY growth length are: λ S = 2 πk ( ω S ) , (35) c S = ω S k ( ω S ) , (36)Λ S = Λ ( ω S ) . (37)Finally, the Lagrangian growth rate is defined as: ς S ≡ c S /λ S . (38) Text S4. Time variations of the characteristic sand flux
The Lagrangian growth rate ς and the propagation velocity c are obtained by fittingexponential and linear functions to the temporal variation of the dune height and po-sition, respectively. However, both are also time-dependent through the characteristicsand flux Q and the slope effect (1 /µ )( u th /u ∗ ) . In the fitting procedure to extract ς and c , we neglect this second contribution, but take into account the temporal variations of Q . Following (26) and (27), the Lagrangian growth rate and propagation velocity areproportional to the characteristic sand flux, such that we can write at each time step: ς t = Q t L ¯ ς, (39) c t = Q t L sat ¯ c, (40)where ¯ ς and ¯ c are non dimensional, assumed to be independent of time. Following theexponential growth and the linear propagation of a bump in the linear regime of the October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 13 instability, its height H and position x at a time step t can be expressed as: H t = H exp (Σ t ς t δt ) , (41) x t = x + Σ t c t δt, (42)which can be rewrittten as:ln (cid:18) H t H (cid:19) = h Q t i L ¯ ς Σ t Q t h Q t i δt = h ς i t ∗ , (43) x t − x = h Q t i L sat ¯ c Σ t Q t h Q t i δt = h c i t ∗ , (44)where h·i denotes the time average. We then see that an effective time t ∗ arises: t ∗ = Σ t Q t h Q t i δt t . (45)It allows the extraction of the time average growth rate and propagation velocity (seeFigure S4), in order to be compared to the theoretical predictions, in which the averagecharacteristic sand flux should be used. Text S5. Distributions of the hydrodynamic coefficients and uncertainties
In order to compare the predictions of the spatial linear stability analysis to the fielddata, we need to find the hydrodynamic coefficients A and B . The field data are dis-tributed around central values (Figure 3 of the main article), due to the intrinsic vari-ability in dune fields. It is important to note that linear stability analysis predicts sucha variability: the initial noise is selectively amplified in a range of wavelengths aroundthe maximum growth rate. Here, we perform an inversion analysis and we determinethe most probable values of A and B , given the observations, using a simple method.We consider values of wavelength, growth length, Lagrangian temporal growth rate andpropagation velocity inside the following intervals, which are representative of the peaks October 5, 2020, 2:44am - 14
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY of the distributions: λ ∈ [100 , , (46)Λ ∈ [90 , , (47) ς ∈ [0 . ,
1] yr − , (48) c ∈ [3 . , .
1] myr − . (49)From a quadruple { λ, Λ , ς, c } , we compute the best tuple {A , B} , taking the optimumprediction of equations (35-38). This method results in a distribution for the hydrody-namic coefficients A and B , determined from the statistical variability of the dune field.It exhibits clear dominant values, slightly different from the average ones due to the dis-tribution shape. By taking the number of independent measurements into account, theerror bars on A and B can be graphically represented by the ellipse of figure S5(a). Thetilt of the ellipse clearly shows that the hydrodynamic coefficients covary, which indicatesa better precision on the ratio B / A than on the individual values of A and B .Average and most probable values also depend on the fixed parameters Q , L sat , µ and u ∗ /u th , whose uncertainties are given in Table S1. Taking those uncertainties into accountresults in multiple different distributions (see for example Figure S5b). Note that all thecorresponding ellipses have roughly the same sizes, fixed by the dispersion of the measuredquantities. The ratio B µ / A is robustly determined as it is in principle related to both λ/L sat and Λ /L sat . In particular, the theory predicts that the ratio Λ /λ ∝ A / B µ shouldbe constant, and independent of any other parameter. Figure S6 shows that the field datanicely obey this property. The theory also predicts that the ratio λc should be proportional October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 15 to h Q t iA . As a consequence, any uncertainty on h Q t i results into an uncertainty on A .Figure S6 shows that the ratio λc/ h Q t i is, as predicted, independent of λ .We therefore have two different sources of uncertainty; one resulting from the dispersionof the measured quantities, and one from the uncertainty on the fixed parameters. Thevalue of L sat is required to calculate λ and Λ, and directly affects the estimate of B / A .Likewise, the value of Q is required to calculate c , and mostly affects A . Finally, theamplitude of the slope effect, (1 /µ )( u th /u ∗ ) directly affects B , as B µ remains the same.The resulting uncertainty on hydrodynamic coefficients results from all contributions,leading to the final result: A = 3 . ± . B = 1 . ± . References
Charru, F., Andreotti, B., & Claudin, P. (2013). Sand ripples and dunes.
Annual Reviewof Fluid Mechanics , , 469–493.Courrech du Pont, S., Narteau, C., & Gao, X. (2014). Two modes for dune orientation. Geology , , 743–746.Dur´an, O., Claudin, P., & Andreotti, B. (2011). On aeolian transport: Grain-scaleinteractions, dynamical mechanisms and scaling laws. Aeolian Research , , 243–270.Fourri`ere, A., Claudin, P., & Andreotti, B. (2010). Bedforms in a turbulent stream:formation of ripples by primary linear instability and of dunes by nonlinear patterncoarsening. Journal of Fluid Mechanics , , 287–328.Iversen, J. D., & Rasmussen, K. R. (1999). The effect of wind speed and bed slope onsand transport. Sedimentology , , 723–731. October 5, 2020, 2:44am - 16
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
Ungar, J. E., & Haff, P. K. (1987). Steady state saltation in air.
Sedimentology , ,289–299. October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 17 h/d . . . C o n e s l o p e (d) Figure S1. (a-c) Evolution of the sand pile when the tube is pulled up. Contour detectionis shown by the blue lines. (d) Measured slope as a function of the non dimensional pile height,where d = 435 µ m is the d of the particular sample we have used. The blue line is the meanvalue, and the shaded area represents one standard deviation above and below it. October 5, 2020, 2:44am - 18
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
Figure S2.
In situ image of representative gypsum grains on the upwind margin; their oblateshape and large size can be seen. The contrast and brightness of grains have been adjusted sotheir shapes are identifiable.
October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 19 .
00 0 .
05 0 .
10 0 .
15 0 .
20 0 .
25 0 .
30 0 . L sat k − . − . − . . . . . T e m p o r a l g r o w t h r a t e L s a t Q σ × − σ T . . . . . . . P r o p aga t i o n v e l o c i t y L s a t Q c = L s a t Q ω / k k T c T .
00 0 .
05 0 .
10 0 .
15 0 .
20 0 .
25 0 .
30 0 . L Q ω = L Q kc − . − . − . . . . . Sp a t i a l g r o w t h r a t e L s a t Λ − × − Λ − − W a v e nu m b e r L s a t k × − ω S k S (a)(b) Figure S3.
Dispersion relations corresponding to (a) the temporal linear stability analysisand (b) the spatial linear stability. Purple lines and dots indicates the most unstable mode. In(b), opaque and transparent curves correspond to the positive and negative branches k + and k − ,relating to upstream and downstream propagating waves, respectively. Here, A = 3 . B = 1 . Q = 38 m yr − , L sat = 2 . µ = 0 . u ∗ /u th = 1 . October 5, 2020, 2:44am - 20
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY . . . . . . . t [yr]0123 C o rr ec t e d t i m e t ∗ [ y r ] y = x t ∗ [yr] − C r e s t p o s i t i o n x − x [ m ] t ∗ [yr]01234 C r e s t h e i g h t l og ( H / H ) × − (a)(b) (c) Figure S4. (a) Corrected time as a function of real time. (b) Position of a bump crest withrespect to time. The dashed line is a linear fit giving c = 5 . − . (c) Height log-ratio of abump crest with respect to time. The dashed line is a linear fit giving ς = 0 .
012 yr − . Errorbarscome from the uncertainty due to the procedure for peak detection. October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 21 . . . . . A . . . B (a)01234567 P D F . . . . . A . . . .
40 0 .
45 0 .
50 0 .
55 0 .
60 0 .
65 0 . B / A . . . . . . √ A + B (c)05101520 P D F .
40 0 .
45 0 .
50 0 .
55 0 .
60 0 .
65 0 . B / A . . . . . . Figure S5.
Distributions of the hydrodynamic parameters resulting from the inversion process.(a) Contour lines of the distribution when the parameters Q , L sat , µ and u ∗ /u th are fixed to theirmean value (Table S1). The markers represent the mean value, and the ellipse the confidenceinterval at 95%. (b) Variation of this ellipse when taking into account the uncertainty of theparameters given in Table S1. Blue and green ellipses are for minimum and maximum plausiblevalues of L sat , respectively. Plain and dashed ellipses are for minimum and maximum plausiblevalues of flux Q . Squares and stars are the mean values of the distributions for minimum andmaximum plausible values of the slope effect (1 /µ )( u th /u ∗ ) . (c) and (d) are the same than (a)and (b), but represented in the space spanned by the modulus √A + B and the ratio B / A . October 5, 2020, 2:44am - 22
GADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
50 100 150 200 250Wavelength, λ [m]02468 Λ / λ (a) DataTemporal analysisSpatial analysis 50 100 150 200 250Wavelength, λ [m] − λ c / h Q t i (b) DataTemporal analysisSpatial analysis Figure S6.
Pertinent non-dimensional numbers in the inversion process. (a) Ratio betweenthe wavelength and the growth length. Following the temporal analysis, Λ /λ ∝ A / B µ . (b) Ratiobetween the wavelength and the propagation length h Q t i /c . Following the temporal analysis, λc/ h Q t i ∝ A . The plain lines represent the predictions of the linear stability analysis for theparameter values given in Table S1, with A = 3 . B = 1 . Table S1.
Table of the parameters used for the linear stability analysis at the White Sandsdune field. The star denotes direct measurements, while the others quantities are derived frommeasurements using calibrated laws. The air density is the value at 20 ◦ C, and the correspond-ing error includes corrections for elevation as well as for temporal variations in humidity andtemperature. Parameter Notation Value UnitGrain diameter* d ± µ mGrain bulk density* ρ p ±
100 kg / m Avalanche slope* µ . ± . ρ f . ± . / m Saturation length L sat . ± . u ∗ /u th . ± . Q ± / yr October 5, 2020, 2:44am
ADAL ET AL.: SPATIAL AND TEMPORAL DEVELOPMENT OF THE DUNE INSTABILITY
X - 23
Table S2.
Measured and predicted values of the characteristics of the emerging dune pattern.The peak value corresponds to the maximum of the smoothed distribution. The peak width isestimated at 3 / λ Growth length Λ Lagrangian growth rate ς Propagation velocity c (m) (m) (yr − ) (myr − ) Data from global method
Peak value 122 112 0.013 4.4Peak width 19 70 0.036 1.5
Data from ‘peak-to-peak’ method
Peak value 116 86 0.015 4.6Peak width 80 125 0.06 3.2
Predictions under a unidirectional wind (red line in Figure 3)
Mean 115 245 0.032 7.6Minimum 105 204 0.018 5.9Maximum 134 330 0.049 10.1
Predictions with correction due to reversing winds (orange line in Figure 3)