An Idealized 3D Model of Interfacial Instability of Aeolian Bedform

: An idealized morphodynamic model is constructed for formation of the aeolian sand ripples from small bottom perturbations of a two-dimensional sand bed. The main goal of the analysis is to evaluate the influence of the gravity flow (including “impact-induced gravity flow” in the reptation flux and “topography-induced gravity flow” in the creep flux) on the formation of the aeolian sand ripples and to clarify the relative contribution of various factors to the bed instability. A 3D linear stability analysis reveals that gravity flow appreciably affects the dynamics behaviors of aeolian sand ripples, which decreases the growth rate of sand ripples, tends to stabilize the sand bed, and leads to longer wavelength. We found that the competition between the destabilizing effect of reptation flow and the stabilizing effects of gravity flow leads to pattern selection. The along-crest diffusion of topography driven by impact and gravity is beneficial to the transverse stability of sand ripples, producing sand ripples with straighter and more continuous crests. For moderate values of D , the most unstable mode has zero value of the transverse wavenumber ( k y = 0), thus corresponding to aeolian ripples with crests perpendicular to the wind. Moreover, when the impact angle is 9–16°, it has little effect on the characteristics of sand ripples for the initial stage of ripple development. For every increase of the impact angle by 1°, the initial wavelength only increases by about 1.5%. In conclusion, the influence of the gravity flow on the dynamics of sand ripples formation stage cannot be neglected.


Introduction
Aeolian sand ripples are dynamic small-scale regular bedforms observed in sandy environments on Earth and other planets [1][2][3]. They are periodic morphological patterns with wavelengths ranging from a few centimeters to dozens of centimeters and heights of less than 1 cm, the profile of sand ripples is asymmetric. In general, aeolian sand ripples migrate, with their crests almost perpendicular to the prevailing wind direction, at a typical rate of about tens of millimeters per minute [4]. Predicting the behavior of aeolian sand ripples has practical interest because they tend to destroy farmland and cause desertification [5].
The process-based models have been developed to analyze the dynamics of aeolian sand ripples. The studies of the processes that lead to the formation of sand ripples has shown that these regular patterns arise as instabilities of the morphodynamic system describing the interaction between topography and sand transport [6][7][8][9][10][11][12].
The model of aeolian ripples can be traced back to the work of Bagnold [13], who proposed the wavelength is determined by the characteristic of saltation length. Sharp [14] found a relation between the wavelength and shear velocity, but the saltation length is larger than the wavelength. The key role of reptation in aeolian ripples formation was suggested by Anderson [6]. In short, the ripples are generated due to slope-dependent impact intensity of the saltating grains that eject reptating grains; the ripple's lee side is then screened from the most energetic impacts, which makes the spatial differences of the reptation flux. Hoyle and Mehta [7] further considered grain rolling or avalanching under gravity. Yizhaq et al. [9] extended the continuous model to 3D structures on the sand surface, refined the reptation flux model, and predicted a saturation of ripples wavelength long term. Recently, based on direct numerical simulations, an instability mechanism leading to the formation of aeolian sand ripples was first suggested by Durán et al. [10]: resonant grain trajectories, whose length is close to the initial wavelength and whose splash causes the mass transport to the ripple crests. However, this study did not involve the later stage of ripples development. Sharp [14] suggested that the wavelength depends on the height of the sand ripples and on the impact angle of the saltation particles. When the impact angle is smaller, the shadow zone length becomes larger, producing longer sand ripples. With wind tunnel experiment and numerical simulation, Schmerler et al. [11] found that the shadow zone mechanism was more dominant in the subsequent nonlinear stage in which successive merging increased aeolian ripples. However, its relationship with the initial wavelength was not clearly demonstrated.
The gravity flow includes impact-induced gravity flow and topography-induced gravity flow, which affect the transverse instability of aeolian megaripples [15]. The impact of the topography-induced gravity flow on the contribution of aeolian bedforms has been reported by Hoyle and Metha [7], while the impact-induced gravity flow reduces the flux upslope and increases the flux downslope [15]. According to wind tunnel experiments and field observations, researchers gradually clarified the influence of wind speed and other factors on the formation and evolution of sand ripples [4,13,14,[16][17][18]. The grain size distributions play an important role in determining the size of the bedform, which has been the subject of considerable investigations [13,14,16,17,[19][20][21][22][23][24][25].
Although significant progress has been made both in the prediction of sand ripples formation and geometric characteristics, much remains to be done. So far, most studies have focused on the one-dimensional sand bed, that is, sand ripples are invariant in the direction transverse to the wind. However, natural sand ripples often show transverse instability [15], resulting in a 3D pattern, which has been rarely investigated [9]. Quite a few previous studies have neglected the effect of transverse sand transport on sand ripples. In addition, the effects of gravity flow on the characteristics of sand ripples are still not well involved in a process-based sand ripples model, and the relative contribution of various factors to the bed instability has not investigated systematically.
In the current study, a 3D idealized model is constructed to systematically investigate the initial stage of sand ripples formation. The model is based on the study of the linear stability of the flat sand bed. The main goal of the analysis is to evaluate influence of the gravity flow on the formation of the aeolian sand ripples and to clarify the relative contribution of various factors to the bed instability. The main innovation is to systematically analyze the above processes, including separately and in combination, then further investigate the transverse instability of the ripple crests and consider the gravity flow on dynamics of sand ripples.

Bed Evolution
The evolution of the sand bed is described by the mass conservation equation for the grains [6,9]: where is the porosity of the sand bed (typically = 0.35) and is the total horizontal volumetric flux per unit width in the x-and y-directions. Equation (1) simply relates erosion/deposition processes to spatial increases/decreases of sand transport.

Sand Transport
Sand transport is the sum of three contributions: saltation, reptation, and creep. The saltating grains only provide the energy for the system and do not directly contribute to the formation of sand ripples [6,9]. We therefore assume that the saltation flux is spatially homogeneous (i.e., ∇ • = 0).
The saltation impact causes the ejection of particles from the sand bed; the ejected particles, which have much lower energy than the saltating particles, make single small hop along the surface of the sand bed. The saltating particles are highly energetic and continue in saltation after rebounding from the sand bed [26].
We assume that the saltating particles impact the sand bed surface at an angle to the x-axis and to the y-axis (see Figure 1). The impacts cause erosion of sand bed, while the saltating grains continue to rebound from the sand bed. The erosion rate ( , , ) of sand grains ejected at ( , ) per unit time is proportional to the impact intensity, which depends on the surface orientation with respect to the direction of saltation, based on the geometrical consideration [19,27], giving  To evaluate the reptation flux per unit width over a horizontal bed, we use the approach of Yizhaq et al. [9]; thus, where ( ) is the probability distribution of the reptation length and is chosen from a normal distribution [19]; the parameter ( = , ) represents the change in the reptation flux due to the slope [9,15]. As for the creep flux of grains caused by gravity-driven, the sand grains tend to roll in the direction of steepest descent under the influence of gravity; for the three-dimensional sand ripples, the gravity-induced velocity of the sand grains can be written as [12] = − sin (4) where g is the acceleration of gravity and r is a function of the grain packing; in the present work, for simplicity, we set r to be constant. The horizontal creep flux of sand grains is proportional to the horizontal speed of rolling, cos ; hence, the horizontal creep flux of grains is where R is the equivalent partial thicknesses of the creep layer.

Scaling Procedure
To simplify the problem but also indicate the possible scenarios, we choose = , = , i.e., the primary transport is in the x direction and the secondary transport is in the y direction, where is the small parameter. The following dimensionless variables are introduced: where is the mean reptation length. After dropping the tilde decoration, Equation (1) becomes:

Linear Stability Analysis
The stability of the planar sand bed is investigated by a normal mode analysis. We introduce small perturbations of the sand bed (strictly infinitesimal), hence the problem can be linearized and the small-amplitude sinusoidal topographic perturbations of the form can be considered: where ℎ ( ) is the perturbation amplitude, is a small quantity, and = + is the topographic wave vector. Substituting Equation (8) into Equation (6), using the linear approximation and neglecting all nonlinear terms, keeping only the linear terms we obtain the perturbed bed evolution equation where ( ) and ̂ are the Fourier transform of ( ), = ( ) . The solution to (9) is obtained by writing Then, we obtain the complex growth rate of the bed perturbation as a function of the topographic wave vector; its real and imaginary parts can be written as: The complex growth rate contains two pieces of information: the real part gives the growth or decay rate of the amplitude of perturbation (i.e., by R( )), and the imaginary part describes the propagation of perturbation in the space (i.e., by Im( )). The bed is linearly stable, if the amplitudes decay for all bed perturbation modes (R( ) < 0), whereas if at least one mode bed of the amplitude grows (R( ) > 0), the bed is unstable.
In order to investigate the mechanisms responsible for the sand ripples dynamics, it is necessary to identify the various contributions to the complex growth rate. The growth rate Equation (11) consists of two contributions involving the creep flux and reptation flux in the x direction and the y direction, which depend on the slope of the bed and on the impact of saltation, respectively. It is the sum of the growth rates due to creep and reptation, according to These contributions are given by In turn, the reptation growth rate Re( ) can be written as This consists of the contributions due to the impact-induced reptation flow Re , , and impact-induced gravity flow Re , . These contributions are given by:

Results
The values of the parameters representative of the typical sand transport conditions on Earth are considered. The thicknesses of the creep layer is related to the grain size and the number of rolling or sliding grains per unit surface, ~ , where is the number of rolling or sliding grains per unit surface. The erosion rate is ~ , where is the number of impacting grains per unit time and surface and is a dimensionless ejection rate coefficient. Following Durán et al. [28], the erosion rate can be expressed as = ⁄ , where is saltation flux and is the mean saltation length. In our model, the reptation length distribution ( ) is extracted from the study of reference [19].
We first analyze the growth rate Re( ), which depends on the creep and is negative and proportional to and ; therefore, the slope-dependence of the creep flux stabilizes short-wavelength perturbations, that is, the effect of the gravity dampens the largest wave numbers, which favors a smoothing of the bedform surface. The growth Re( ) contains two contributions, Re , and Re , ; Re , is positive, proportional to cot and cot , and related to the probability distribution of the reptation length. The smaller the angle of impact, the larger cot ( = 1,2) becomes; it plays a destabilizing role. Whereas Re , is negative and stabilizes the sand bed, it is proportional to ( = , ); the parameter represents the change in the reptation flux due to the slope effect. Physically, the saltating grain affects the sand bed and induces some grains to roll or slide down the slope. This mechanism is termed an "impact-induced gravity flow" and acts to reduce the flux upslope and to increase it downslope [15]. Hence, the pattern selection results from a compromise between the physical mechanisms of stabilizing and destabilizing. The competition between these two opposite mechanisms results in a region of amplified in the space , . By changing the relative strength of these two effects, we can control the instability of the system and approximate a marginal stable state. Figure 2 shows the growth rate Re( ) of the bed perturbations as a function of the topographic wave numbers and , in which (A) the total growth rate is the sum of the contributions of (B) reptation and (C) creep. In turn, the contributions to the growth rate due to reptation are given by the sum of (D) impact-induced reptating grains flow effect, which is the growing mechanism, and (E) impact-induced gravity flow effect, which is the decaying mechanism. The results show that the reptation is the growing mechanism and the slope effect is the decaying mechanism.
The numerator represents the intensity of creep transport on the inclined plane, and the denominator is the horizontal flux of reptation grains on the flat bed, so the ratio between the creep flux and the reptation flux is ∝ , in which measures the ratio of the creep flux to the reptation flux, which responsible for the relative transport intensity of the two transport modes. The next important issue is to determine the range of possible values for the model parameter D in real physical situations. Note that the experimental determination of reptation and creep is rather inadequate, especially the creep transport driven by gravity. It is therefore instructive to study the influence of these two transport modes on the marginal stability curve. In Figure 3, we find that the unstable region is very sensitive to a change of the constitutive parameter. The results show that the formation and evolution of the sand ripples occur under the combined action of saltation, reptation, and creep; the saltating grains energize the reptating population; and the gravity causes sand grain creep. So, it is proved that the instability of the sand bed is determined by the competition of the reptation and creep. For small values of , such as < 0.01, the effect of the impact-induced gravity flow is much stronger than the topography-gravity flow (i.e., ( = , ) ≫ ); the topography-induced gravity flow is less than 1% of the flux of reptating grains, and it has almost no effect on the stability of sand ripples. Therefore, the shortwavelength stabilization originates from impact-induced gravity flow. In contrast, for large values of , the topography-induced gravity flow plays a leading stabilizing role; a band of unstable modes appears only if reptation flow prevails over topography-induced gravity flow (or more precisely, if is greater than a critical value given by = , see Figure 4). In order to efficiently investigate the effect of the parameter D on the aeolian ripples formation, a one-dimensional investigation has been made by fixing = 0. The marginal stability conditions shown in Figure 4 were plotted in the space ( , ). The critical value Dcr of D, above which sand ripples do not appear, and sand bed stable can be easily estimated from the graph, which also shows that the sand bed becomes unstable in a finite range of wave numbers kx once D is less than its critical value. The maximum growth rate of sand bed perturbation versus its wavenumber kx for different values of D is represented by the red line shown in Figure 4. The maximum wavenumber tends to an infinitesimal value ( → 0) as D tends to Dcr, and the wavelength of preferred mode tends to infinity ( → ∞), so the sand bed is stable, and the sand ripples do not appear. For = 0, the creep flux disappears; the model is similar to the model of Yizhaq et al. [9].  for different values of D. We find that the characteristics of these contour plots highly depend on the value of . The zero growth contours on the graph locate the marginal instability curve. For = 0.1, the perturbation component is characterized by the maximum growth rate, i.e., the component that dominates the subsequent development of sand bed instability is characterized by , = (2.9, ±2.3) (see Figure 5B). These wavenumbers correspond to the wavelength of about 2 cm. For = 1, compared with the last case, it is clear that the instability region (i.e., regions with positive growth rates) decreases (see Figure 5C); it reflects the fact that the most unstable growth rate reaches 0.64 at , = (1.5, 0) and the initial wavelength is about 4 cm. The initial wavelength scale is similar to the observed wavelengths [4] and simulates wavelengths that are about 25 times the grain diameter [10]. It is also found that the most unstable mode is longitudinal (i.e., ≈ 0) and therefore corresponds to aeolian ripples having crests perpendicular to the wind. When = 1.5, the instability region shrinks to a small contour, almost becoming a single point (see Figure 5D). Hence, with the increase of , the marginal instability curve asymptotically approaches a point; if continues to increase, the sand bed becomes stable and sand ripples will not be formed. We investigate the effect of the impact angle for the sand ripples formation. The saltating grains typically impact the bed at angles between 9° and 16° [29][30][31]. In doing so, we fix the topographic wave vector in y direction as 0, that is, = 0. Figure 6 shows the growth rate Re( ) of the sand bed perturbations versus for different values of versus of . When the impact angle is increased, the growth rate Re( ) shows that the flat sand bed is more stable. Moreover, the value of , which gives rise to the maximum value of Re( ), increases, thus showing that bigger impact angle generates longer sand ripples; the wavelength of the preferred mode increases slightly with the increase of impact angle at the initial stage (see Figure 6 (insets)). This possible reason is that the sand ripples are very small at the initial stage of sand ripple development; the greater the impact angle, the more effective the saltation impact, which will produce greater the reptation flux so that the sand ripples will grow. However, for the subsequent nonlinear growth stage, the shadow effect can be one of major mechanisms contributing to the growth of sand ripples; the smaller the impact angle, the larger the shadow zone and the longer the wavelength [11,14]. The influence of impact angle on the dynamics of sand ripples is different in the initial linear stage and the later nonlinear stage. It remains for future further studies to determine the influence of impact angle in different stages of ripples formation. It is interesting to discuss the unstable transverse mode, although the mode does not prevail in the linear regime; since its growth rate is smaller than that of the preferred mode, we may wonder whether they play an important role in the transverse development of sand ripples. The influence of topography-induced gravity transport on the transverse instability is illustrated in Figure 7A; it shows that for small values of (for example = 0.01), a wavy crestline with finite wavelength is formed, which corresponds to the most unstable mode . This is the early stage of sand ripple formation, and the ripples are small (i.e., D is small), so the lateral transport of gravity flow is small, which reduces the lateral coupling and produces sand ripples with more irregular and lateral instability. This is consistent with the characteristics of the early formation stage of sand ripple observed in wind tunnel and field. With the increase of D, most unstable transverse mode is → 0. For example, when = 0.5, we get ≈ 0. It reflects that the most unstable mode is a longitudinal mode in both modes; thus, the corresponding sand ripples has more straight and continuous crests perpendicular to the wind direction. In fact, as the sand ripples grow, the gravity flow increases (i.e., increases). The lateral transport of gravity flow has the capability of reducing along-crest topographic irregularities and producing sand ripples with straighter and more continuous crests perpendicular to the wind. This is also consistent with the features observed in wind tunnels and field. Figure  7B plots the growth rate versus for = for different values of ; it shows similar results, where measures the strength of the sand flux by impact-driven along the y direction, which is responsible for the lateral coupling [15]. A larger lateral (along-crest) sand transport tends to stabilize the transverse mode, and the wavelength of the transverse mode becomes larger, which has the capability of reducing along-crest topographic instability and irregularities and producing bedform-like structures with straighter and more continuous crests. The numerical simulations of Yizhaq et al. [15] support our analysis results that owing to the large lateral flux, the along-crest height of sand ripples is almost the same; the transverse tends to be stable, reducing along-crest ripples sinuosity, creating a straighter crestline. In fact, the aeolian sand ripples are generally straight and continuous, which requires the physical coupling process along the crest of the sand ripples, and that lateral (along-crest) sand transport provides effective physical mechanisms to couple the topography at different along-crest locations and to overcome local instability. In the proposed model, the lateral transport is mainly contributed to by the along-crest diffusion of topography driven by impact and gravity; with the increase of and , the lateral transport increases, which provides a mechanism for extending sand ripples in the along-crest direction. The normal ripples composed of monodisperse fine sand particles have larger transverse transport, which can fill downwind concavities and flatten downwind convexities; the sand ripple is straight and continuous [32]. However, for megaripples, coarse particles accumulate at the crest, which reduces the lateral sand transport along the ripple crest and increases the irregularities and transverse instability along-crest [15]. The transverse aeolian dunes are not straight, while the origin of transverse instability is due to different drift velocities caused by the different heights along the cross-wind direction [33] and due to the fact that the lateral transport on the dune slip face prevails over that on the upwind face [34]. This difference may arise from the strong scale separation between the grains' average hop length and the typical size of the dune, while there is no such scale separation for smaller-scale sand ripples. The transverse modes do not grow as fast as the longitudinal modes in the linear regime, but they could couple each other in the nonlinear regime. A further nonlinear analysis should be made to investigate the interaction between different modes and predict the defect behavior and bedform straightening in the non-transverse cases.

Discussion
Aeolian sand ripple is a typical geomorphic type resulting from the interaction between topography and sediment transport, which is mainly caused by the impact on a sand bed of saltating grains entrained by the wind. In this work, we propose an ideal model to systematically investigate influence of the gravity flow on the formation of the aeolian sand ripples, to discuss the lateral stability of sand ripples and to clarify the relative contribution of various factors to the bed instability; this is rarely seen in the existing theoretical model. Based on grain-scale simulations, Durán et al. [10] argued that the instability is driven by resonant grain trajectories with a distance close to the initial ripple wavelength, and the stability effect originates from the slope-dependent sand transport. However, it was not clearly explained how the stabilizing is related to the phase-locking caused by particles jumping from one crest to another, and the lateral stability of sand ripple is not involved. In this work, we discussed the stability mechanism through the diffusion transport of gravity flow and analyzed the influence of gravity flow on the formation of sand ripples. We found that the lateral gravity flow provides effective physical mechanisms to overcome local transverse instability and to produce straight crests. However, due to the assumption that the space of saltation flux is uniform in continuous model, the model has some limitations. The present 3-D model is based on the linear stability stable analysis; it cannot be used to investigate the saturated amplitude attained by the growing perturbation for the later stage. Indeed, with the growth of sand ripples and reach large values, nonlinear effects become strong and the linear analysis fails. It remains for future studies to investigate the nonlinear regime.

Conclusions
We have developed a new idealized 3D aeolian sand ripples model that is able to systematically investigate gravity flow and other factors related to processes in aeolian ripples dynamics in the formation stage. We found that gravity flow (including "impactinduced gravity flow" and "topography-induced gravity flow") decreases the growth rate of sand ripples significantly. The impact-driven reptation flow has a destabilizing effect, while the gravity flow tends to stabilize the sand bed and leads to longer wavelength. On the other hand, when the impact angle is increased from 9° to 15°, the initial wavelength is increased by only 0.2 cm, which means the impact angle has little effect on the characteristics of sand ripples at the initial stage of ripple development. The lateral gravity flow provides effective physical mechanisms to overcome local transverse instability and to produce straight crests. Finally, the results show that the relative transport intensity of the two transport modes essentially effects the instability of sand bed and consequently affects the formation of sand ripples.
The model provides a simplified description of the actual phenomenon in the sand ripples formation; however, it takes into account the main processes that affect the formation of sand ripples (e.g., saltation, reptation, and creep). This ideal model is, of course, far from natural situations; however, we expect to gain an improved understanding of the instability mechanisms of the three-dimensional bedforms. Indeed, this can be considered as a canonical situation. In the future, some aspects need to be further improved. For example, a better model capable of describing the mesoscale structure of aeolian sand transport could be introduced to study the effects of wind speed on the behavior of aeolian sand ripples, to overcome the shortcomings of the existing continuous model that cannot describe the dependence of the wavelength on wind speed.