_________
* Corresponding author e-mail address: takurosuzuki@ffpri.affrc.go.jp
Application of an MPS-based model to the process of debris-flow
deposition on alluvial fans
Takuro Suzuki
a,*, Norifumi Hotta
b, Haruka Tsunetaka
a, Yuichi Sakai
baForestry and Forest Products Research Institute, Japan, 1 Matsunosato, Tsukuba-shi, Ibaraki 3058687, Japan bThe University of Tokyo, 1-1-1 Yayoi, Bunkyo-ku, Tokyo 1138657, Japan
Abstract
A modified moving particles simulation model (MPS-DF) to simulate inundation and sediment deposition of debris flows is presented. This model is based on the moving particles semi-implicit (MPS) method, which was originally used for incompressible viscous fluid flows with free surfaces. In the MPS-DF model, the constitutive equations of Egashira is introduced to the MPS method. In Egashira’s theory, debris flows are treated as a continuum and sand grains are expressed using sediment concentration. Thus, each particle has a variable sediment concentration value. In this study, we tested the applicability of the MPS-DF model for the formation process of alluvial fans. For this purpose, flume experiment was conducted. The experimental flume consisted of a straight channel 6.0m long and 0.1m wide, with an inclination of 15°, connected to an outflow plain. The inclination of the outflow plain decreased gradually from 12° to 3°. At the straight channel, 5.0m long erodible bed with a thickness of 0.2m was present. Water was supplied from upper end for 60 s. at the rate of 3,000 cm3/s and debris flow was generated by entraining the erodible bed. Debris flow inundated and deposited sediment at the outflow plain and an alluvial fan was formed. Numerical simulations were also performed with the MPS-DF as well as a depth-integrated method based on the shallow water equations (2D simulation). 2D Simulation results of alluvial fan shape and flooding area were laterally spread and significantly different from those of experiment. The results of the MPS-DF were more similar to experimental results. Natural channels and lateral levees were formed as well as experiment. However, the alluvial fan shape of MPS-DF was slightly wider than that created during the experiment. This is thought to be due to the behavior of pore water of deposied layer, such as the seepage of water out of the deposited layer once the deposition process has been completed.
Keywords: numerical simulation; particles method; alluvial fan formation process
1. Introduction
Debris flows cause enormous damage. To mitigate the damage, it is important to predict their range of influence. Numerical simulation is an effective tool for predicting debris-flow inundation (e.g., Liu et al., 2012; Pastor et al., 2014). In Japan, residential areas are sometimes built on alluvial fans at valley exits of steep mountain rivers. An example that illustrates this situation is a debris-flow disaster that occurred in Hiroshima in August 2014 (Nakatani et al., 2017). In this disaster, it was suggested that houses and other structures that existed in the alluvial fan area at the exit of the valley impeded the flood and deposit of debris flow. It is also shown that when debris flow occurs continuously in multiple streams, the fan formed by preceding debris flows affects the fan formation process of subsequent debris flows (Chen et al., 2016). In numerical simulations based on shallow water equations, it is necessary to average flow velocity distribution and sediment concentration distribution in the vertical direction, and it is therefore considered difficult to reproduce such complicated behaviors (Suzuki and Hotta., 2015, 2016).
In recent years, there have been advances made in research on particle simulation methods for debris flow based on Smoothed Particle Hydrodynamics (SPH) (Monaghan, 1988) or Moving Particle Simulation method (MPS) (Koshizuka and Oka, 1996). Most of these studies treat debris flow as one viscous fluid (e.g., Laigle et al., 2007; Wang et al., 2016). With these methods, separation of sand grains and water is not reproduced. Suzuki and Hotta (2015,
2016) developed a particles simulation method for debris flow (hereinafter referred to as MPS-DF) that introduced the constitutive law for flow resistance and sediment concentration of debris flow (Egashira et al., 1989, 1997) to the MPS method. With this method, it has been shown that the erosion and deposition process can be reproduced in constant width flumes (Suzuki and Hotta, 2015, 2016). Conversely, the process of debris-flow deposition on alluvial fans involves transverse spreading. In this study, we tested the applicability of the MPS-DF model for the deposition process on alluvial fans. For this purpose, flume experiments were carried out using gravel of uniform grain size and numerical simulations were also performed. The numerical simulations were performed using the MPS-DF method as well as a method based on the shallow water equation, and the results of both methods were then compared.
2. Experiment
2.1.Experimental flume and materials
The experimental flume consisted of a 6.0 m long and 0.1 m wide straight channel, with an inclination of 15 ̊,
connected to an outflow plain (Fig. 1). The gradient of the outflow plain changes by 3 degrees every 1 m, and the 3-degree area only has a length of 2 m. The transverse direction is horizontal. In the straight channel part, sand grains with a depth of 20 cm were deposited in the 5 m downstream part. The height of the sediment surface was equal to the height of the connecting outflow plain. The average grain diameter was 0.265 cm while the specific gravity was 2.6.
2.2.Method and measurement
Water was supplied at 3,000 cm3/s from the upstream end for 60 sec, while debris flow was generated through the
erosion of sediment in the straight channel. The sediment concentration was about 32% at the beginning and about 16% at the time of 60 sec, when it was measured using the erosion depth roughly estimated from the side-view videos of the flow. The debris flow separated into water and sand grains at the outflow plain, forming an alluvial fan. The formation process was photographed using three digital cameras from the top of the outflow plain. The photographed images were analyzed using SfM software (Photoscan Professional, Agisoft LLC) to create point clouds. Elevation data was then created by processing the point clouds and known coordinates, while the deposition depth was calculated as the difference from the original altitude data (De Haas et al., 2014).
3. Numerical Simulation
3.1.Outline of MPS-DF method
Suzuki and Hotta (2015, 2016) developed the MPS-DF method by introducing the constitutive equations developed by Egashira et al. (1989, 1997) into the MPS method (Koshizuka and Oka, 1996). The outline is as follows.
In the MPS-DF method, in order to treat the debris flow as a one-fluid model, the debris flow itself is divided into aggregates of particles that possess parameters of sediment concentration. First, Suzuki and Hotta (2015, 2016) modified the resistance law of debris flow as follows, so that it can be introduced into the framework of the MPS method. = + ( ) + 2 ( ) (1) = ∗ ( ) − ∗ ! " # $ %&− 1 ( )* + ,-. /0 (2) ( ) = 123+ 2 4 5, 23= 89(1 − :5) %$& !, 2 = 8 ( )! ! (3)
where τc represents the shear stress, τcy represents the yield stress, ρw =1.0 is the density of water, c represents the
sediment concentration, u represents the flow velocity, z represents the axis in the depth direction perpendicular to the flow direction of the debris flow, Kd is the coefficient of particle collision, Kf is the coefficient of pore water turbulence,,
c* is the sediment concentration in the bed, σ =2.65 is the gravel density, /s = 38.5° is the friction angle, d is the
diameter of the sediment particles, e = 0.775 is the coefficient of restitution, and where kg = 0.0828 and kf = 0.16 are
empirical constants (Egashira et al., 1989, 1997). The first term of Eq. (1) is the yield stress, and the second and third terms show the dynamic stress for each differential order. θ is an angle formed by the flow velocity vector of the particle and the horizontal vector. τc that was obtained by converting τc into a vector in the rectangular coordinate
system was introduced instead of the viscous term used in the MPS method.
Next, a model was constructed in which the value of the sediment concentration moves among neighboring particles. First, the equilibrium concentration gradient, gce is calculated by substituting the above parameters of a focused
particle (its number is defined as i) into the concentration distribution formula developed by Egashira et al. (1989, 1997). Then, a concentration gradient with neighboring particles, gc is calculated. Fig. 2 shows the concept of the variable sediment concentration model. The change in the sediment concentration of a focused particle was calculated for the underlying particles, because the concentration distribution was obtained by integrating the concentration gradient from the riverbed to the water surface. The particle shown in Fig. 2 tends to result in an increase in c because
gce > gc. Therefore, ccp(i), which expresses the magnitude and direction of the change in c, is defined as follows:
;(<) = =>(( ?− ( ) (4)
Here, l0 is the standard particle distance. The time-derivative of c is derived using the kernel function, w(r),
assuming that c changes in proportion to the difference in ccp(i) between neighboring particles. Here, w(r) is defined by Koshizuka and Oka (1996).
@ A=B CD∑ F(G)AIJ 5( ;(<) − ;(H)) (5)
Here, n0 is the standard particle number density and is a constant parameter of the MPS method (Koshizuka and
Oka, 1996), r is the inter-particle distance, and T is the relaxation time. Notably, T needs to be small enough to satisfy the local equilibrium of sediment concentration (Suzuki and Hotta, 2015, 2016). Therefore, we adopted the method proposed by Suzuki et al. (2016) to link the relaxation time with the increment time.
With the MPS-DF method, for example, when debris flow moves from a steep slope to a gradual slope, the sediment concentration for the upper layer particles moves to the lower layer particles. This is because the equilibrium concentration, calculated using a low gradient, is lower than the current value. The sediment concentration moves to reduce the difference between the current high sediment concentration and equilibrium concentration. When the sediment concentration of the lower layer particles exceeds a certain level, the yield stress term represented by Eq. (2) becomes larger than the shearing force, and the particles are in an immovable state; that is, they become deposited particles (a state in which vibration is slightly repeated). This process corresponds to a deposition process, while the opposite process corresponds to an erosion process.
Fig. 2. A schematic diagram illustrating the concept of the variable sediment concentration model (Suzuki and Hotta, 2016). 3.2.Treatment of momentum in the concentration change model
Moving the sediment concentration is equivalent to exchanging gravel with water. Because gravel and water have different specific gravities, moving the sediment concentration alone cannot satisfy the momentum conservation law. In previous studies, the influence on the calculation result was small, so we ignored the non-conservation of momentum. However, there is a possibility that it may affect the planar deposition process. Therefore, when calculating the change in sediment concentration using Eq. (5), the change in momentum corresponding to a movement in sediment concentration was calculated simultaneously (hereafter it is called the modified MPS-DF). In this study, numerical simulations were performed using a method that does not consider the change in the momentum (called MPS-DF in this study) and the influence of momentum conservation was verified.
3.3.Two-dimensional simulation based on shallow water equations
For comparison, a two-dimensional simulation based on the standard shallow water flow equations (Hereafter, it is called the two-dimensional simulation) was also performed. This method adopted Egashira's constitutive law for flow resistance (Egashira et al., 1989, 1997). Eq. (1) - (3) are obtained by differentiating the resistance law used in this method with z. This method also adopted erosion rate formula as follows (Miyamoto and Ito, 2002).
K = ,-.(+ − +?)|M| (6)
,-.+?= ($ %($ %⁄⁄&& )O) ,-./0 (7)
Here, E is the erosion rate, and θe is the equilibrium gradient corresponding to the sediment concentration.
According to Eq. (6) and Eq. (7), erosion and deposition are determined only by the relationship between the equilibrium gradient corresponding to the average sediment concentration in the vertical direction and the bed gradient.
3.4.Calculation condition
In the MPS-DF method, the particle diameter was 0.5 cm. In the two-dimensional simulation, an orthogonal grid with a width of 5 cm was used. These values were determined from the relationship between resolution and calculation time. The other parameters mentioned above are material properties or empirical constants.
4. Results
4.1.Alluvial fan formation process
The experimental results and calculation results of the alluvial fan formation process are shown in Fig. 3, while the final deposition depth is shown in Fig. 4. In Fig. 3 and 4, the outlet of the straight channel was set as the 0 m point.
Fig. 3. Alluvial fan formation process. The yellow broken line indicates a region with no surface flow.
The time at which the front of the debris flow arrives at the outlet of the straight channel was set to 0 sec, and the results are shown every 15 sec.
In the experiment, the deposit expanded linearly in the downslope direction for approximately 3 m without spreading horizontally. After that, the flow direction changed to the right side. Then, when the deposition on the right side progressed to a certain extent, the flow changed to the left side direction. The final result was the formation of an almost symmetrical alluvial fan. A waterway was formed near the outlet of the straight waterway. A waterway was formed near the outlet of the straight waterway, similar to the result obtained by De Haas et al, in which a waterway and a natural levee were formed near the outlet of the straight waterway (De Haas et al., 2015).
In the result obtained from the two-dimensional simulation, deposition began immediately from the outlet of the straight channel that gradually spread in the downstream side and transverse direction. Compared to the experimental results, the spread in the transverse direction was large and the distance covered in the downward flow direction was short. It can be said that the deposition process/shape from the two-dimensional simulation was significantly different from that obtained from the experimental results.
In the result obtained using the MPS-DF method, a linear deposition shape was formed in the initial process, but it expanded in the transverse direction from around the point beyond the 1 m point. The flow direction did not change and a symmetrical fan was formed. Similar to the experimental results, a waterway and natural levee were also formed. The final fan shape spread in the transverse direction, in contrast to the trend obtained from the experimental results.
In the result of the modified MPS-DF, the lateral spread was suppressed to some extent due to the movement of the momentum. However, as time passed, it expanded in the lateral direction, and its deviation from the experimental results increased.
4.2.Final deposition depth
Fig. 4 shows the result of the final deposition depth and Fig. 5 shows the longitudinal section of the final deposition depth at the center of the outflow plain. The experimental result shows only the region analyzed using the SfM software. Comparing the results of the MPS-DF and modified MPS-DF models, it was observed that the result of MPS-DF model showed a slight downstream deposition. Regarding the point with the largest deposition depth, the results obtained from the two-dimensional simulation were significantly different from the experimental results. Conversely, the MPS-DF and modified MPS-DF results generally agreed with the experimental results. In particular, the result of the modified MPS-DF was nearly quantitatively consistent with the experimental results.
5. Discussions and Conclusions
In the two-dimensional simulation, deposition occurred immediately near the outlet of the straight channel and spread laterally. In the experiment, deposition gradually progressed linearly from the straight channel outlet. Only the water in the upper layer of the flow spread sideways. In the two-dimensional simulation, the movement of sand grains between meshes was calculated using the sediment discharge, obtained by multiplying the flow rate by the transport concentration (Miyamoto and Ito, 2002). The flux sediment concentration which is obtained by dividing sediment discharge by total discharge is a mesh-specific value; that is, it is the same value in the downward flow direction and the lateral direction in the two-dimensional simulation. Therefore, it is believed that the situation in which the upper layer of water selectively spreads sideways was not reproduced, and instead, that the deposition widely spread in the lateral direction. In the MPS-DF and modified MPS-DF methods, since the direction of movement and sediment concentration of each particle was calculated, it is possible to automatically evaluate sediment transport in the downward flow and lateral directions. Therefore, the evaluation was greatly improved from the two-dimensional simulation.
Furthermore, in the two-dimensional simulation, the deposition distance in the downward flow direction was short. Additionally, in the two-dimensional simulation, an erosion/deposition rate equation was used, but the momentum change process was not considered. Therefore, the fact that the deposition occurs suddenly corresponds to the sudden loss of momentum by the sand grains. The absurdity of momentum accompanying erosion/deposition in the calculation using shallow water flow equations was pointed out by Iverson and Ouyang (2015), and is a difficult problem to solve. Since the motion equation was solved for each particle in the MPS-DF method, the deviation from the experimental results was small. However, the non-conservation of momentum occurs when the momentum transfer accompanying concentration movement is ignored. In the case of the deposition process, the concentration moved from the upper layer with a large momentum to the lower layer with low momentum, implying that the momentum was gradually lost. As a result, this may have allowed for the easier deposition of sediment. The modified MPS-DF that solved this problem yielded results that more closely matched the experimental results, especially in the early stage of deposition. In the early stage, the deposition was formed in the upper stream in the MPS-DF model, more so than in the modified MPS-DF model, while final deposition was formed in the most downstream side. This is considered to be the effect of sediment re-erosion by subsequent debris flow with low concentration, because the eroded sediment becomes easy to move for the reason opposite that of the deposition process. However, the longitudinal deposition results of the modified MPS-DF model were almost in agreement with the experimental results.
As described above, the modified MPS-DF model has the highest reproducibility for the longitudinal deposition results. To accurately reproduce the erosion/deposition process, it is important to strictly evaluate the sediment concentration distribution, flow velocity distribution, and momentum conservation law. The particle method is effective for conducting this evaluation. However, in the later process, the results obtained using the modified MPS-DF model deviated from the experimental results in terms of lateral spreading. One of the reasons for this was that bedload was generated from the edge of the deposit and continued to flow downstream. As the size of the calculation particles was larger than the average particle size of the sand grains, it proved impossible to represent the movement of individual sediment particles. Although reducing the size of the calculation particles could solve this problem, it is not realistic in terms of the calculation load. Therefore, other solutions need to be considered.
Another reason for the deviation of the results is due to the behavior of seepage flow in the deposited sediment. Fig. 3 shows the region where surface flow did not occur. The vicinity of the surface in this region is unsaturated. In other words, water moves differently than sand grains. Even in the saturated region, seepage flow seemed to occur in the deposited sediment. Since water selectively flowed out from the pore, it was believed that sand grains were difficult to spread laterally. It is difficult to reproduce the behavior of seepage flow using the particles method based on the one-fluid model. However, it is possible to reproduce unsaturated seepage flow by giving parameters of water content ratio to particles and moving them among neighboring particles. Based on this model, it is necessary to improve the calculation model so as to reproduce the seepage of water out of the deposited layer.
Acknowledgements
References
Chen, H.X., Zhang, L.M., Gao, L., Yuan, Q., Lu, T., Xiang, B., Zhuang, W.L., 2016, Simulation of interactions among multiple debris flows: Landslides, v. 14, p. 595-615.
De Haas, Ventra, D., Carbonneau, P. E., Kleinhans, M. G., 2014, Debris-flow dominance of alluvial fans masked by runoff reworking and weathering: Geomorphology, v. 217, p. 165–181.
De Haas, T., Braat, L., Leuven, J.R.F.W., Lokhorst, I.R., Kleinhans, M.G., 2015, Effects of debris flow composition on runout, depositional mechanisms, and deposit morphology in laboratory experiments: Journal of Geophysical Research: Earth Surface, v. 120, no. 9, p. 1949-1972. Egashira, S., Ashida K., Yajima H. and Takahama J., 1989, Constitutive equations of debris flow: Disaster Prevention Research Institute Annuals,
Kyoto University, no. 32 B-2, p. 487-501.
Egashira, S., Miyamoto, K., Ito, H., 1997, Constitutive equations of debris flow and their applicability: Proc. 1st International Conference on Debris-Flow Hazards Mitigation, ASCE, p. 340-349.
Iverson, R. M. and Ouyang, C., 2015, Entrainment of bed material by Earth-surface mass flows: review and reformulation of depth-integrated theory, Reviews of Geophysics, v. 53, p. 27-58.
Koshizuka, S. and Oka, Y., 1996, Moving Particle Semi-implicit Method for Fragmentation of Incompressible Fluid: Nuclear Science and Engineering, v. 123, p. 421-434.
Laigle, D., Lachamp, P., Naaim, M., 2007, SPH-based numerical investigation of mudflow and other complex fluid flow interactions with structures: Computational Geoscience, v. 11, p. 297–306.
Liu, J., Nakatani, K., Mizuyama, T., 2012, Hazard Mitigation Planning for Debris flow Based on Numerical Simulation Using Kanako Simulator: Journal of Mountain Science, v. 9, p. 529–537.
Miyamoto K. and Ito T., 2002, Numerical simulation method of debris flow introducing the erosion rate equation: Journal of the Japan Society of Erosion Control Engineering, v. 55, no. 2, p.24-35.
Monaghan J. J., 1988, An introduction to SPH: Computer Physics Communications, Volume 48, Issue 1, p. 89-96.
Nakatani, N., Kosugi, M., Satofuka, Y., Mizuyama, T., 2017, Influence of housing and roads on debris flow flooding and deposition in alluvial fan areas : Case study on debris flows in Hiroshima, Japan, in August 2014: Journal of the Japan Society of Erosion Control Engineering, v. 69, no. 5, p. 3-10.
Pastor, M., Blanc, T., Haddad, B., Petrone, S., Sanchez Morles, M., Drempetic, V., Issler, D., Crosta, G. B., Cascini, L., Sorbino, G., Cuomo, S., 2014, Application of a SPH depth-integrated model to landslide run-out analysis: Landslides, v. 11, p. 793–812.
Suzuki T. and Hotta N., 2015, Development of modified particles method for simulations of debris flows based on constitutive equations and its application to deposition process: Journal of the Japan Society of Erosion Control Engineering, v. 68, no. 1, p. 13-24
Suzuki, T. and Hotta, N., 2016, Development of Modified Particles Method for Simulation of Debris Flow Using Constitutive Equations: International Journal of Erosion Control Engineering, v. 9, no. 4, p. 165-173.
Suzuki, T., Hotta, N., Iwata, T., 2016, Accelerating methods for particle model of debris flows and its application to entrainment process: Proceedings of the 8th Symposium on Sediment-Related Disasters, p.193-198.
Wang, W., Chen, G., Han, Z., Zhou, S., Zhang, H., Jing, P., 2016, 3D numerical simulation of debris-flow motion using SPH method incorporating non-Newtonian fluid behaviour: Natural Hazards, v. 81, p. 1981–1998.
