This paper presents a new meshless numerical scheme to overcome the problem of shock waves and to apply boundary conditions in cases of dam-break flows in channels with constant and variable widths. The numerical program solves shallow water equations based on the discrete mixed subdomain least squares (DMSLS) meshless method with collocation points. The DMSLS meshless method is based on the minimization of a least squares functional defined as the weighted summation of the squared residuals of the governing equations over the entire domain and requiring the summation of residual function to be zero at collocation points in boundary subdomains. The collocated discrete subdomain meshless method is applied on the boundary, whereas the collocated discrete least squares meshless technique is applied to the interior domain. The meshless scheme extends for dam-break formulation of shallow water equations. The model is verified by comparing computed results with analytical and experimental data for constant and varying width channels. The developed model is also used to study one-dimensional dam-break problems involving different flow situations by considering changes to the channel width, a bumpy channel with various downstream boundary conditions, and the effects of bed friction and bed slope as source terms on wave propagation. The accuracy of the results is acceptable.
INTRODUCTION
Dam-break flows propagate along rivers and can lead to devastating floods, damage to property, and loss of human life. Dam-break problems have been the subject of scientific research by many hydraulic scientists and engineers. Many efforts have been made to understand the mechanism of dam-break flows by conducting analytical, experimental, and numerical studies. Whereas dam break is associated with large free-surface deformation and fragmentation with the interfacial flow and mixing process, it makes modeling of the flow structure very complicated (Biscarini et al. 2010; Bellos et al. 2012; Ozmen-Cagatay & Kocaman 2014).
The most important feature of dam-break equations is that they allow discontinuities and smooth solutions. These discontinuities, often called ‘shock’ or ‘shock waves’, initiate the failure of a number of classical numerical methods. The results of dam-break modeling show that shallow water equations are logically appropriate for capturing discontinuous free surface, and methods for shallow water problems fit sufficiently with experimental and analytical results (Alcrudo & Garcia-Navarro 1993). Shallow water equations are typically used to represent the hydrodynamics and behavior of dam-break problems.
Earlier studies have been primarily based on analytical solutions of shallow water equations for idealized conditions. For example, Stoker (1957) developed an analytical solution to predict dam-break flows in an idealized channel in which the bed slope was assumed to be zero and the friction term was ignored. Chanson (2006) proposed an analytical solution of dam-break waves with flow resistance, and then it was applied to simulate tsunami surges on dry coastal plains. However, analytical solutions for shallow water equations which are nonlinear hyperbolic partial differential equations (PDEs) cannot be obtained except in special and simple cases. With the development of computer technology, hydrodynamic models based on one-dimensional (1-D) and two-dimensional (2-D) shallow water equations are increasingly being used in numerical solutions for dam-break flows. Numerous studies in the literature solved these equations using different numerical techniques, such as finite difference, finite volume, and finite element methods based on the mesh system (Bellos et al. 1991; Stansby et al. 1998; Mohamed et al. 2003; Janosi et al. 2004; Bukreev & Gusev 2005; Liang et al. 2006; Ghadimi & Reisinezhad 2012; Duricic et al. 2013; Di Francesco et al. 2015). However, mesh-based numerical methods have difficulty dealing with complicated problems such as in dam-break flows with violent free surface and interfaces. The main difficulty in using such methods is the capturing of sharp interfaces or shock waves. Furthermore, the complicated procedure needs to specify those cells containing the free surface. Also, in mesh-based methods, the problem of numerical diffusion is common due to advection terms (Zhihua 1989).
Recently, the hottest interest in the field of computational fluid mechanics has been concentrated on development of the next generation of numerical meshless methods, which are advantageous in applications with discontinuous free surface or large interfacial deformations compared to mesh-based methods. As their name implies, one common characteristic of all these methods is that they do not require the traditional mesh to construct the numerical formulation; they require node generation instead of mesh generation. In other words, there is no pre-specified connectivity or relationship among the nodes; thus, the computational costs associated with mesh generation are greatly reduced. The computational advantages of a meshless method suggest that they have the potential to solve a broad class of scientific and engineering problems (Liu 2002; Liu & Gu 2005). Over the last decade, meshless methods have increasingly been used for the solution of PDEs to solve both shallow water equations and Navier–Stokes equations. Ferrari et al. (2010), Kao & Chang (2012), and Ataie-Ashtiani et al. (2008) used the smooth particle hydrodynamic (SPH) concept to solve both shallow water equations and Navier–Stokes equations. The moving particle semi-implicit (MPS) method has been applied for simulation of fixed-bed dam-break flow on dry-bed and wet-bed downstream channels (Shakibaeinia & Jin 2010; Khayyer & Gotoh 2010). Darbani et al. (2011) used the natural element method to simulate 2-D shallow water equations in the presence of strong gradients. Despite their advantages, meshless methods have the difficulties of imposing boundary and treatment of boundary values, whereas the unreliable selection of weight functions and the complexity of algorithms for computing interpolation functions are all serious technical problems in such methods. This is in contrast to mesh-based methods which often have this property. Consequently, the imposition of boundary conditions requires certain attention in mesh-free methods and may degrade the convergence of the method (Firoozjaee & Afshar 2010a). Special techniques, such as the penalty method by Zhu & Atluri (1998), the transformation of approximate nodal values to actual nodal values by Cai & Zhu (2004), the nodal integration method by Beissel & Belytschko (1996), and the efficient computation of shape functions by Beitkopf et al. (2000) have been proposed to overcome these problems.
Recently, a family of collocation-based meshless methods has been emerging in the literature. The advantage of collocation methods is their efficiency in constructing the final system of equations because integration is not required, and shape functions are evaluated at nodal positions only. However, the accuracy and robustness of collocation approaches are weak points, especially if the approximation is based on a set of randomly scattered points. Collocation methods may lead to large numerical errors in these cases and involve numerical stability problems. An interesting discussion of these aspects may be found in Fries & Matthies (2004), where certain conditions are proposed which are often not met by standard collocation mesh-less methods. Furthermore, boundary conditions are an issue for collocation methods. The definition and location of the boundary surface may not be easy, e.g., for free surface flow, and methods of applying boundary conditions are not always straightforward. Firoozjaee & Afshar (2010b) proposed the collocated discrete least squares meshless (CDLSM) method to solve steady-state incompressible Navier–Stokes equations. Shobeyri & Afshar (2009) discussed simulating free surface problems using the discrete least squares meshless method. Afshar et al. (2009) used the CDLSM method to solve transient and steady-state problems.
In CDLSM, the penalty method is used and can only approximately satisfy the boundary conditions. Accuracy is affected by selection of the penalty coefficient; these boundary conditions are only approximately imposed, depending on the magnitude of the penalty coefficients. Theoretically, the larger the penalty coefficients, the more accurate the enforcement of boundary conditions will be. It is difficult to choose a set of penalty factors that are universally applicable for all kinds of problems. One hopes to use the largest possible penalty factors, but penalty factors that are too large often give numerical problems, and trials may be needed to choose a proper penalty factor (Liu 2002; Liu & Gu 2005).
In this study, the authors propose a new method named collocated discrete subdomain meshes (CDSM) method on the boundary for more accuracy as explained in the section ‘Implementation of boundary conditions’. In this approach, the residual functional can be written in summation form and minimized at collocation points; the sum of the residuals must be zero for each boundary subdomain. It is noteworthy that the boundary subdomains are usually centered by nodal points on the boundary.
This study investigated the discrete mixed subdomain least squares (DMSLS) meshless method using collocation points. To generalize this method, an attempt was made to check the balance of this approach for some dam-break flow simulations (where the definition and location of the boundary surface is not easy) for different bed and geometry conditions with maximum variables in boundary conditions (where applying boundary conditions are not straightforward). To verify the prediction results, this technique was tested for different experimental dam-break flows considering bed friction, bed slope, bumpy channels with various downstream conditions and varying width channels, and well-known analytical solutions for rectangular channels of constant width. Finally, the results were summarized and conclusions were drawn. The results of this study illustrate the potential ability of the DMSLS method to reproduce the detailed features of flow structure for highly transient flow with the large geometrical changes of the computational domain.
GOVERNING EQUATIONS
NUMERICAL SCHEME
Local support domains in the meshless method.
Moving least squares method
Approximate function and nodal parameters
in the MLS approximation.
. The coefficients can be obtained by minimizing the following weighted discrete norm:, and is the nodal parameter of u at . Equation (10) is a functional, weighted residual that is constructed using the approximated values and the nodal parameters of the unknown field function. Because the number of nodes used in the MLS approximation is usually much larger than the number of unknown coefficients , the approximated function does not pass through nodal values. In this research, the cubic spline weight function explained in more details by Liu (2002) and Liu & Gu (2005) was considered.
MLS shape functions generally do not satisfy the Kronecker delta condition. Hence, cannot be treated like nodal values of the unknown function
, the dimension of the support domain is determined by is the dimensionless size of the support domain, and is the nodal spacing near the point at . Generally, leads to good results for many of the problems that have been studied (Liu 2002; Liu & Gu 2005). It should be noted that the support domain is usually centered by a point of interest. It should also be noted that the polynomial basis of the second order was used in this study for 1-D problems to construct MLS shape functions, such as 1-D dam-break cases. Therefore, the accuracy of this approach is the second order in space.Mathematical formulations
can be obtained through interpolation: is the number of nodal points having spatial coordinates in their support domain.
can be obtained through interpolation:
Implementation of boundary conditions
is residual at collocation points inside the Dirichlet subdomain. Each Dirichlet subdomain is usually centered by a nodal point on the Dirichlet boundary, and is the number of collocation points inside each Dirichlet subdomain.
defined as:.
NUMERICAL CASES
is the celerity is the Courant number that should be less than 1.
norm is used: and are the simulated water depth and the exact solution at the ith point, respectively.Idealized dam-break problem in rectangular channel
, the dam is broken instantaneously. A positive wave propagating downstream and a depression wave traveling upstream can be observed. The initial upstream water depth of the dam is 10 m, and the corresponding initial downstream water depths are given as 0.0 m (dry bed) and 2.0 m (wet bed) for each case. The analytical solution for dry bed on slope channel is given as (Chanson 2006): is the initial water depth of the dam break, H is the water depth, is the distance from the dam break, and t is the time. More details about the analytical solution on a wet, friction bed can be found in the literature (Stoker 1957; Chanson 2006).
Numerical solutions of (a) water depths and (b) velocities for frictionless, dry bed.
Numerical solutions of (a) water depths and (b) velocities for frictionless, dry bed.
Numerical solutions of (a) water depths and (b) velocities for frictionless, non-dry bed.
Numerical solutions of (a) water depths and (b) velocities for frictionless, non-dry bed.
as shown in Figure 5. The comparison between numerical and analytical solutions shows an acceptable agreement.Numerical solutions of water depths for dam break on sloping bed.
, height of , width of , a bottom slope of 0.93%, initial length of water column , and initial heights of water column upstream , respectively. Table 1 summarizes the norms between DMSLS and MPS with respect to experimental data. As can be seen, the results of this new approach (DMSLS) are in good agreement with the experimental data and have greater accuracy compared with the MPS method.The norms between the simulated water depths and exact solutions
| Reference | DMSLS | MPS |
|---|---|---|
| norm | 0.0488 | 0.0584 |
| Reference | DMSLS | MPS |
|---|---|---|
| norm | 0.0488 | 0.0584 |
Comparison of DMSLS meshless method with MPS and experimental data (Ying et al. 2004).
norms between the numerical and exact solutions of different NNP and NCP, and results demonstrate that the norm decreases as NNP increases. As shown in Table 2, collocation points instead of nodal points can be increased to obtain better results.The norms between the simulated and exact water depths for different numbers of nodal points (NNP) and numbers of collocation points (NCP) and water depth ratios
| Water depth ratio | ||||||
|---|---|---|---|---|---|---|
| | |||||
| Case | No. of nodal points | No. of collocation points | CPU time (s) | L2 | CPU time (s) | L2 |
| A | 101 | 201 | 6.62 | 0.083 | 7.26 | 0.306 |
| B | 401 | 8.36 | 0.052 | 9.06 | 0.053 | |
| C | 1,001 | 13.59 | 0.033 | 15.05 | 0.025 | |
| D | 201 | 401 | 59.42 | 0.049 | 64.71 | 0.051 |
| E | 1,001 | 78.09 | 0.031 | 84.77 | 0.032 | |
| F | 2,001 | 114.06 | 0.022 | 122.59 | 0.024 | |
| Water depth ratio | ||||||
|---|---|---|---|---|---|---|
| | |||||
| Case | No. of nodal points | No. of collocation points | CPU time (s) | L2 | CPU time (s) | L2 |
| A | 101 | 201 | 6.62 | 0.083 | 7.26 | 0.306 |
| B | 401 | 8.36 | 0.052 | 9.06 | 0.053 | |
| C | 1,001 | 13.59 | 0.033 | 15.05 | 0.025 | |
| D | 201 | 401 | 59.42 | 0.049 | 64.71 | 0.051 |
| E | 1,001 | 78.09 | 0.031 | 84.77 | 0.032 | |
| F | 2,001 | 114.06 | 0.022 | 122.59 | 0.024 | |
Numerical solutions of water depths for various numbers of collocation points under water depth ratio of 0.2.
Numerical solutions of water depths for various numbers of collocation points under water depth ratio of 0.2.
It is notable that increasing nodal points to obtain better results yields an increase in the dimension of the coefficient matrix, but increasing the number of collocation points yields better results without increasing the dimension of the coefficient matrix. Thus, this method can be easily used in practical cases with high efficiency and less computational effort, and it can simulate dam-break problems in real size (in the field) rather than on a small scale with good accuracy and less simulation time.
Real dam-break problem in a rectangular channel
Numerical solutions of water depths for frictional, dry bed.
Comparison of idealized and real numerical solutions of water depths for non-dry bed.
Comparison of idealized and real numerical solutions of water depths for non-dry bed.
Straight channel with a triangular bump in bed
Straight channel with a triangular bump in the bed. Side elevation showing gauge positions indicated with a suffix number after G which means distance from the gate, i.e., G4 indicates that it is 4 m away from the gate (Morris 2000).
Straight channel with a triangular bump in the bed. Side elevation showing gauge positions indicated with a suffix number after G which means distance from the gate, i.e., G4 indicates that it is 4 m away from the gate (Morris 2000).
Straight channel with a triangular bump in the bed, wet bed with open end. Water surface profile at s.
Straight channel with a triangular bump in the bed, wet bed with open end. Water surface profile at s.
Straight channel with a triangular bump in the bed, wet bed with closed end. Water surface profile at s.
Straight channel with a triangular bump in the bed, wet bed with closed end. Water surface profile at s.
Straight channel with a triangular bump in the bed, dry bed with an open end. Comparisons of experimental and numerical water depths at selected gauges.
Straight channel with a triangular bump in the bed, dry bed with an open end. Comparisons of experimental and numerical water depths at selected gauges.
Straight channel with a triangular bump in the bed, wet bed with an open end. Comparisons of experimental and numerical water depths at selected gauges.
Straight channel with a triangular bump in the bed, wet bed with an open end. Comparisons of experimental and numerical water depths at selected gauges.
Straight channel with a triangular bump in the bed, wet bed with a closed end. Comparisons of experimental and numerical water depths at selected gauges.
Straight channel with a triangular bump in the bed, wet bed with a closed end. Comparisons of experimental and numerical water depths at selected gauges.
Varying width dam-break problems
Comparison of the DMSLS solution with experimental data for : (a) x = 0.0 m, (b) x = 4.5 m, (c) x = 8.5 m, (d) x = 13.5 m, (e) x = 18.5 m.
Comparison of the DMSLS solution with experimental data for : (a) x = 0.0 m, (b) x = 4.5 m, (c) x = 8.5 m, (d) x = 13.5 m, (e) x = 18.5 m.
Comparison of the DMSLS solution with experimental data for : (a) x = 0.0 m, (b) x = 4.5 m, (c) x = 8.5 m, (d) x = 13.5 m, (e) x = 18.5 m.
Comparison of the DMSLS solution with experimental data for : (a) x = 0.0 m, (b) x = 4.5 m, (c) x = 8.5 m, (d) x = 13.5 m, (e) x = 18.5 m.
CONCLUSIONS
The moving least squares (meshless) scheme has been developed to solve the 1-D shallow water equations for the computation of dam-break induced flows. A key feature of the model is the use of collocated points to minimize the sum of the squared residual functional. The problem domain is discretized using some nodes. Then, using the MLS shape function with a least squares technique, an asymmetric stiffness matrix is constructed. Sampling (collocated) points are defined to minimize the sum of the squared residual functional over interior domain and to be zero residual function at their boundaries, and field nodes (nodal points) are defined to construct MLS shape functions. Considering the number of sampling (collocated) points to be more than the number of field nodes will improve the results and guarantee stability.
Nonlinear 1-D dam-break problems have been solved using the newly proposed CDSM method for boundary conditions and the CDLSM method for interior domain, and the results have been presented and compared to the analytical and experimental results. The variety of dam-break problems involving different flow regimes by changing the width of the channel, considering the effects of bed friction and bed slope as source terms on wave propagation, and bumpy channel with reflective boundary were modeled. The results clearly indicate that the proposed DMSLS method with collocation points for boundary conditions is capable of obtaining stable and more accurate results in dam-break problems, especially for wet/dry interfaces, different bed geometry conditions, and shock and hydraulic break waves.
The formulation presented in this paper offers a number of advantages over many existing meshless methods used to solve PDEs. In CDLSM, the penalty method has been used for boundary conditions; these boundary conditions were imposed only approximately, depending on the magnitude of the penalty coefficients. It is difficult to choose a set of penalty factors that are universally applicable for all kinds of problems. One hopes to use the largest possible penalty factors, but penalty factors that are too large often give numerical problems, and trials may be needed to choose a proper penalty factor. Thus, this paper introduces the CDSM method to enforce boundary conditions in moving least square mesh-free methods. Another advantage of this method is the use of collocated points. It is notable that increasing nodal points for obtaining better results yields to an increase in the dimension of the stiffness matrix, but increasing the number of collocation points yields better results without increasing the dimension of the stiffness matrix. This issue prevents a large stiffness matrix, increases the efficiency of the scheme, and consumes less time to solve the system of equations. These advantages can help researchers simulate flood events in the field induced from prototype dam break. This method can also be used for a moving bed such as dam-break modeling over an erodible bed, an issue currently under research.
%20cropped.png?versionId=7142)