Modeling and Numerical Simulation of Sand Dunes Formation

We are interested this paper in the modeling and numerical simulation of sand dune formation (on the effect of wind). We propose a simple model which will be approached by a classical numerical method in order to verify validity of our model to describe the phenomena of dune formation. Some numerical solutions of dune evolution and approximation error table are presented to verify that the simplified model is adopted to study the formation of sand dunes.


Introduction
Dunes are accumulations of sands of larger or smaller sizes and constituted of relatively mobile grains on the surface. The dunes appear anywhere sand is exposed to a strenght (wind, water, etc.). The sand dune formation is due to two main phenomena : the effect of fluid on the sand which causes the sand grains to move and the generation of the avalanche phenomena which arises when the slope of the dune exceeds certain critical angle varying between 20 and 40 (Daniel, 2001).
In the last century many works so experimental as theoretical has been carried out on the formation of sand dunes .One of the first reference work in the field is certainly the famous book of Bagnold which dates back from 1941 (Bagnold, 1941). Since then, a great effort of measurement and modeling has been done (Schwammle & al, 2004;Kroy & al, 2002;Schwammle & al, 2005;Hersen, 2004;Matteo & al 2011;Jakson and Hunt, 1975) . Indeed numerical modeling, based on a better understanding of the formation of sand dunes can help in the practical development of effective fighting means. On the other hand, it is a phenomena that remains difficult and complex to model.
In this paper we are interested in modeling and simulating the evolution of a sand dune in time. More clearly, we propose a simplified model from existing models (Andreotti, Claudin, & Daoudy, 1975) that we solve by classical numerical methods.The objective of this study is to show that the simplified model can be used as a model allowing to account for the dynamics of the formation of the dunes and also to emphasize especially the different parameters determining in the analysis of the process of formation of sand dunes.
The outline of this paper is as follow: In the second section we present the general model from which we elaborate our simplified model which will be presented just after, in the third section we proceed to the approximation of our model and we present in the last section some numerical simulations. We end this work with a conclusion and perspectives.

Classical Model
Let Ω be a bounded domain of R occupied by the dune for an instant T > 0. The general mathematical model from which we derive our simplified model is governed by the following equations (Matteo, Luca, & Fabio, 2011): Where: • h (t, x) denotes the height of the dune at every point x of the space and at each instant t; • q (t, x) represents the flux of sand grains transported at any point and at any time; • q sat (t,x) is saturation flux ; • Q sat (t,x) represents saturation flux on flat bed (h=0); • γ denotes the angle of stability of the dune beyond which the avalanche phenomenon occurs; • u * denotes the shear stress velocity; • U * is the shear velocity on a flat bed; • C B : constant parameter; • A and B are phenomenological parameters.
Equation (2.1) models the conservation of mass which describes the interaction between wind and sand sediment.We impose that the slope of sand surface cannot exceed the angle of repose during evolution else avalanche occurs.
The last equation (2-4) models the flux of transportd sand. It is a charge equation that connects the reel flux of sand to the satured sand flux (Schwammle & Hermann, 2005).

Simplified Model
The resolution of the system (2.1) -(2.4) is done by generally using first the Fourier transform technique to approximate equation (2.2). For this reason we develop a simplified model that not only circumvents this difficulty but is more numerically exploitable and must be able to take into account the important aspects of modeling. We then reformulate the system (2.1) -(2.4) in the following way Where α and β are phenomenological parameters which reflect the dependence of the flux on the curvature and the slope respectively and D is characteristic length of dune.

Approximation of the Model (2.5) -(2.10)
To approximate the model, we consider the spatial domain Ω = [a, b].
Let N x be a non-zero natural integer, we define the discretization step of the domain Ω We also subdivide the time domain [0, T] into N t subintervals [t j , t j+1 ] with N t a non-zero integer and t j = j × ∆t for j ∈ {0, 1, · · · , N t }.
We obtain the following system defined at the point ( We approach the partial derivatives by the following formulas.We use a forward scheme of order 1 to evaluate the temporal derivation and a second order centered scheme for the second derivative in space. By substituting these expressions in the system (3.1) -(3.7) we obtain the following explicit Euler scheme for i = 1, · · · , N x − 1; j = 1, · · · , N t − 1 Let us denote by H j , the vector whose components are the approximate solutions at the nodes x i and at the time t j .Then we can rewrite schemes (3.10) in the following matrix form: whose components are the conditions at the edges of the domain.
By approaching the derivative at x by the formula leads to the following implicit Euler scheme for the approximation of our model: We can rewrite scheme (3.13) in the following matrix form : A convex combination of the explicit and implicit Euler schemes with a coefficient t = 1/2 gives us the Crank-Nicholson scheme for the approximation of our model given by the following relations: Rewrite in the matrix form the system become:  The above figures describe the evolution of a sand dune with different approximation schemes. For the starting form of the dune, we chose a polynomial bump of degree 2, h 0 = 2(1 − x 2 ) with a maximum height of 2 m. Figures(1)-(6) illustrate the evolution of a sand dune for different values of discretization of the spatial and temporal domain. Through the different simulations we observe an increase in the height of the dune as a function of the temporal and spatial discretization parameters. It should also be noted that for the implicit and Crank-Nicolson schemes the results are better when the number of discretization points increases while the explicit scheme produces a solution which deteriorates as the number of points increases.The explicit Euler schemes gives solution that diverge when certain conditions are note satisfied.This divergence is shown through figure (5A) which clearly show that the solution provide by the explicit method does not agree with the reference solution and those obtained by the two other methods when the discretization are ∆t = 0.001, ∆x = 0.05. This scheme imposes a restriction on the discretization steps in time and space which is well known by the Courant, Friedrichs, Lewy condition (CFL). We can conclude that the Crank-Nicolson and implicit schemas are better suited for the approximation of our model problem than the explicit method.
Figures (7) and (8) give a simulation of the height of the sand dune as a function of the phenomena parameters α and β.
The Tables 1−5 illustrate the approximation errors obtained by the different schemes for different numbers of discretization points at a time T of evolution of the dune and a given value of the ration ∆t ∆(x) 2 . It is constant that on all the tables, for all the schemes, for a fixed number of points, the error decreases as the iteration time increases. On the other hand, an analysis of Tables 1 and 3 shows that over the same period of evolution of the dune, and for the same ratio ∆t ∆(x) 2 , when the number of points of iteration is increased, the approximation becomes better. When the ratio ∆t ∆(x) 2 is small in front of a condition CFL which is ∆t ∆(x) 2 < 1 2αD , the explicit Euler scheme gives a better approximation than the other schemes. On the other hand, when the ratio is greater than the CFL condition, the explicit Euler scheme explodes.The positive exponents in tables 4 and 5 simply show that the error generated by the explicit Euler method tends to infinity which expresses the nonconvergence of the method and the explosion phenomenon of the given solution by explicit Euler method when the CFL condition is not satisfied. The implicit scheme and Crank-Nicholson always gives better approximations which is perfectly mathematical agreement of the stability performed above.Indeed, it has been shown that the implicit Euler and Crank schemes are unconditionally stable while the explicit Euler scheme is stable if the ratio ∆t ∆(x) 2 is lower than the CFL condition.

Conclusion
In this paper we have presented a simplified model of sand dune formation that we have approached by classical numerical methods and we simulated this model. These simulations allowed us not only to choose the appropriate method for the analysis of this model but above all to numerically determine the determining parameter in the formation of sand dunes. This parameter will enable us in a future work to control the movement of these dunes.