Numerical solution by meshless method of a fully-coupled bed load and shallow water flows

The one-dimensional shallow water equations (SWEs), coupled with bed load (Exner equation) were solved in this paper for the numerical simulation of the sediment transport problem. By using an efficient meshless method based on radial basis function (RBF) that approximates the unknown variables on a set of collocation nodes and Runge–Kutta fourth-order scheme to approximate the temporal derivative. Besides, the artificial viscosity to avoid oscillations, we implement the proposed numerical model of bedload and set up a series of numerical tests using various nontrivial solutions to confirm the efficiency of the suggested approaches. The numerical simulations demonstrate the accuracy, computational efficiency, and ability of the suggested numerical methods in solving the coupled Exner shallow water problem. . * Corresponding author: marwan.19.bens@gmail.com ICEGC'2021 E3S Web of Conferences 336, 00003 (2022) https://doi.org/10.1051/e3sconf/202233600003 © The Authors, published by EDP Sciences. This is an open access article distributed under the terms of the Creative Commons Attribution License 4.0 (http://creativecommons.org/licenses/by/4.0/).


Introduction
In the last decades the sediment transport has become a major topic of interest in the hydraulics community because this phenomenon (i.e. sediment transport) is of great economic importance, for example: to predict the possible bed form changes of rivers and estuaries; to predict the risks of scouring of bridges, weirs and channel banks; to estimate the siltation of a reservoir upstream of a dam wall. In the most numerical simulation of sediment transport problems, the mathematical model is composed of a hydrodynamical component that describes flow continuity and momentum conservation, coupled with a morphodynamical component that expresses sediment continuity. The Saint-Venant equations, often known as shallow water equations, are used to predict the hydrodynamic behaviour of shallow water flows, while the Exner equation is used to compute sediment transport and simulate bed morphology changes. To model the solid transport sediment flux, different equations can be found in the literature such as the Grass equation [1], Meyer-Peter and Müller equation [2], and Van Rijn's equation [3], etc. these models can generally be obtained by empirical methods. The governing equations obtained are formulating a system of nonlinear partial differential equations (PDEs), which in general have no analytical solutions, those equations have usually been solved using classical numerical methods based on mesh generation, such as finite difference, finite volume, and finite element methods, but in recent decades the finitevolume method has been paid significant attention to model the bed load flows. Delis et al. proposed a numerical method based on finite volume relaxation approximation for computing the bed-load sediment transport in shallow water flows, in one and two space dimensions. They modelled the water flow by the wellknown nonlinear shallow water equations which are coupled with a bed updating equation, using a secondorder MUSCL-TVD method for the advection stage while an implicit-explicit Runge-Kutta scheme solves the relaxation stage [4]. There is also Benkhaldoun et al. whose discuss the application of a finite volume method to morphodynamic models formulated for the quasi-steady approach and the coupled approach on unstructured triangular meshes. they used a modified Roe's scheme incorporating the sign of the Jacobian matrix in the morphodynamic system with a wellbalanced discretization was used for the treatment of source terms [5]. A new proposed numerical method combines the method of characteristics with a finite volume discretization of the system. The combined method is simple to implement and accurately resolves the governing equations without relying on Riemann problem solvers, the mathematical formulation of these models consists of the shallow water equations for the hydraulics, an advection equation for the transport of sediment species, and an Exner equation for the bedload transport. The coupled problem forms a onedimensional hyperbolic system of conservation laws with geometric source terms [6]. Bilanceri et al. consider a physical model based on the shallow-water equations coupled with the Exner equation to simulate the morphodynamic flows using implicit time advancing combined with two finite-volume methods [7]. Coupling between the shallow water and the Exner equations has also been rather extensively studied in the literature (see for example [8] and [9]).
In our work, a numerical model based on meshless methods is applied for solving the sediment transport equations. Which is a new class of numerical methods that are free of mesh and they are widely used in recent years to solve numerical problems (see for example [10], [11], and [12]). One among those meshless methods is the radial basis function (RBF) method, extended initially by E.J. Kansa [13]. This method will be presented in our work for the simulation of a onedimensional sediment transport problem. For this purpose, the shallow water equations along with the Exner equations are solved using the RBF method, Besides, the artificial viscosity to avoid oscillations, while the Runge-Kutta fourth-order scheme is used for time discretization. We organize our paper as follows: In section 2 we present the governing equations for the sediment transport problems. In Section 3 the numerical formulations of the meshless method are presented. In Section 4 we present the numerical approximation of the Exner shallow water model by using the RBF method. In section 5 we use some numerical examples in order to check the accuracy and the efficiency of the proposed method, some conclusions are summarized in this section as well.

Mathematical formulation
The shallow water equations (SWEs) are derived from the Navier-Stokes equations under the following simplifying assumptions. First, it is assumed that vertical accelerations are sufficiently small to be neglected and that the fluid is homogeneous; the velocity profile can therefore be considered uniform through the depth. For a fixed bed Z = Z(x) and neglecting viscous and Coriolis effects, the onedimensional shallow water equations, containing geometrical source terms due to the bottom and friction topography, can be written as: where h = h(x,t) is the water depth, u = u(x,t) is the depth-averaged velocity component and g is the acceleration due to gravity. For this study, the righthand side of the momentum equation comprises the bed slope. Note that also furthers terms can be considered such as the bed friction and surface wind effect. In this study we suppose that the bed depends on time variable, that is Z=Z(x,t), where the bed load updating is given by the Exner equation: Where p is the sediment porosity which is assumed to be constant. and Qs represent the bed-load sediment transport flux in the x-direction. In the literature, they are many models that can be used to evaluate this flux.
Here, the Grass model was adopted and the flux is written as: A is a dimensional constant (s/m), that encompasses the effects of grain size and kinematic viscosity and is usually determined from experimental data (0 ≤ A ≤ 1), with m being chosen so that 1 ≤ m ≤ 4. The parameter A is related to the strength of the interaction between the sediment transport and the water flow. When A is near zero, a weak interaction occurs between bed-load and water flow. A strong interaction occurs in the case of A is near one. Finally, when A is equal to zero, the bed has no movement and considered solid. For this paper, we take m=3 (see [8]) which gives: In a coupled form, the system of equations can be written as follow: The above system can be writing as a conservative hyperbolic system: 3 Numerical formulations of the meshless method In purpose of formulating our method to solve the governing equations (5) distinct nodes. At X, we can approximate a function such as U(x) in scattered data interpolation as a linear combination of N radial basis functions as follow: 1 , , Where ߣ are the unknown coefficients to be calculated, ‫ݎ‬ = ฮ‫ݔ‬ − ‫ݔ‬ ฮ is the distance between the points ‫ݔ‬ and ‫ݔ‬ , and ߮(ฮ‫ݔ‬ − ‫ݔ‬ ฮ) is the radial basis function. Many radial basis functions exist in the literature (see Table 1) and in the current study we consider the infinitely smooth multiquadric radial basis function defined as: Where 0 H z is the shape parameter controlling the fitting of a smooth surface to the data. The selection of optimal values for the shape parameter has been the subject of several studies, see [14][15][16] and [11] among others. In the present work we used the following selection: Where m d denotes the smallest nodal distance, and N is the number of nodes.
Similarly, for any linear differential operator L 1 , , The expansion coefficients ߣ ‫)ݐ(‬ are obtained by solving the following linear algebraic system of N × N equations of ܷ , that we construct from the approximation (8).
)/ U Where 1 2 , So, the unknown expansion coefficients Λ are calculated by: The spatial derivatives at ‫ݔ‬ are calculated as: Using (18), the matrix forms of (19) can be formulated as: And 1 xx xx xx U D ) ) U U (21) In which ‫ܦ‬ ௫ and ‫ܦ‬ ௫௫ are called the first-order and second-order similarly differentiation matrix. where ߔ ௫ and ߔ ௫௫ represent the matrix of first and second derivatives of radial basis functions.

Using RBF method for numerical approximation of Exner Shallow Water model
Now we will use the method described above in our work for the simulation of one-dimensional sediment transport problem. For this purpose, the shallow water equations along with the Exner equations are solved using the RBF method, Besides, the artificial viscosity in order to avoid oscillations (see for example [19] and [20]), while the Runge-Kutta fourth-order scheme is used for time discretization.
In [21] and [8], several reformulations of the system (5) were considered. By using the product rule

Artificial Viscosity
For particular hydraulic problems involving shock waves, the numerical model is required to represent a steady or unsteady discontinuity, where the presence of oscillations in the solution is expected [19], and it sometimes may grow over time. However, by introducing a small amount of artificial diffusion, it is possible to damp these oscillations [22], [23]. We, therefore, augment the right-hand side of the hyperbolic system (22) with the artificial diffusion term: Where D P is a tunable viscosity coefficient and calculated by its relation to the Peclet number which is defined as follows: It often controls the stability of the numerical solutions.
In the case of convection-dominated flow, the Peclet number is large but finite [24], and the effect of the diffusion term (23) becomes negligible. Therefore, a natural and simple way to stabilize the solution is to reduce the Peclet number.
Using the RBF interpolation of an arbitrary function (21), the augmented term (23)

Applying the RBF method for spatial discretization
So, we applied our meshfree method described above on the right-hand side of the hyperbolic system (22)

Time discretization
To achieve a higher order of accuracy, we use the Runge-Kutta fourth-order scheme. The procedure to advance the solution from the time ‫ݐ‬ to the next time ‫ݐ‬ ାଵ is carried out as: and ‫ݐ∆‬ is the time step. To achieve stability, for this explicit scheme, the time step must meet the following criterion: Where ‫ܮܨܥ‬ is the courant number, such that 0 < ‫ܮܨܥ‬ < 1, and ݀ denotes the smallest nodal distance between collocation points.

Accuracy test
In order to check the accuracy of our numerical method, we solve the problem of lake at rest used to verify the balance between the convective and the bed fluxes knowen as C-property. In this situation, we also purpose that the bed can be moved, and is defined by: > @ The channel is of length 10 m and the initial water height and discharge are defined as: As in [8], [21] the parameters used in this problem are the same used in the one-dimensional sediment transport flux in (4). The numerical solution is computed with N refers to the total number. Where ܻ (ே௨) is the numerical solution, and ܻ (ோ) is the reference solution, and N is the total number of points.  As it is shown in the table every increase in the number of points causes a decrease in errors, and in Fig.1. we notice that there is a good accuracy shown by the numerical results obtained by the RBF meshfree method with a very small error shown in Table 2. So, the RBF method performs well for this problem and produces accurate solutions without requiring a complicated mesh and with a simple space discretization.

Transport of parabolic sediment layer
To study the convergence of the proposed meshless method, we consider the governing equations in (22). The parameters of this study are those used in [8], [21]. Thus, the channel is of length 1000 m and the initial bed is defined as 2 300 sin , 300 500 (0, ) 200 0 , The initial water level and initial velocity are given as: As in [8,21], the parameter m = 3, the porosity ߳ = 0.4, and ܳ is a constant discharge taken the values 10 m 2 /s and 50 m 2 /s. Two values A=0.001 and A=1 are used in the one-dimensional sediment transport flux in (22). The space domain is discretized in 400 nodes for our computations. We have used a reference solution computed in [8] and used in [25] for comparison reasons. The initial conditions obtained can be seen in Fig.  2.
(a) (b) Fig. 2. The initial conditions for the bed and water level.

Weak interaction A=0.001
In this case, the bed is interacting slowly with the water flow. To simulate this, we use the values A=0.001 and Qc=50 m 2 /s at t = 2500 s.
In Fig. 3. we present the results obtained of the water height and the bed at t=2500s, using the RBF meshfree method, the results are compared with reference results computed in [8]. As we can see there is no notable difference between our results and those from [8], which prove that the RBF meshfree method performs very well for this sediment transport problem.

Strong interaction A=1
In this case, the bed is interacting strongly with the water flow. To simulate this, we use the values A=1 and Qc=10 m 2 /s at t = 700 s.
The Fig. 4. presents the obtained numerical results by using the RBF meshfree method of the bed profile and the water height with A=1 at time t=700s and QC=10 m 2 /s, compared with analytical results. We notice that there is a good agreement between the analytical solution and the numerical results obtained by the RBF meshfree method.

Conclusion
In the present study, we have presented an investigation on the MQ-RBFs meshless method for handling sediment transport problems by solving the coupled system of shallow water equations and sediment transport equation. A radial basis function has been used to approximate the space derivative on global collocation points, while the Range-Kutta fourth-order scheme has been used to approximate the temporal derivative of the differential equations. The MQ-RBFs meshless method is an uncomplicated method that simple to implement with good flexibility. It is a meshfree method what makes it an interesting method compared to the classical methods. For the verification of the considered method, we used several standard sediment transport problems with different bed interactions and flow regimes. The RBF mesh-free method performs very well and a good agreement was obtained between the numerical results and the analytical solutions for the transport problem of the parabolic sediment layer, which makes it promising to apply also to other situations and problems.