Unified particle method with application to the 2013 Huangping landslides and tsunamis event of China

. In this study, a unified particle method is presented to simulate the fluid-solid coupling problem in broad range of scales such as landslide and tsunamis. First, a general overview of the method is addressed, and the governing equations are solved in a Lagrangian form. Second, the method is used to simulate the Scott Russell wave generator experiment containing simple fluid-solid coupling, and the relationship between the simulation results and the experimental results is analyzed to verify the validity of the model. Then, the method is applied to the 2D processes simulation of Huangping landslide generated impulse waves in large scale. The results show that: in the verification test, the wave amplitude error between the simulation and the experimental data is almost 0. In the application case of Huangping landslides and tsunamis, the maximum impulse wave height obtained by the simulation is close to the local observations, which indicates that the method has high accuracy and credibility.


Introduction
Impulse waves may be generated by landslides, shore instabilities and rock falls in reservoir area. If the impulse waves is strong enough, it may cause water to pass over the crest or even destroy the building. This will block the river, threaten the navigation of ships, and threaten the safety of life and property of coastal residents [1] . In particular, such geological disasters are more likely to occur in the reservoir area where water storage or water storage is being carried out. On the one hand, the construction of the reservoir will inevitably change the original natural geological conditions of the bank slope of the reservoir area, such as the bank slope structure, and then affect the stability of the slope [2] . On the other hand, the reservoir water will infiltrate and support the unstable geologic bodies and landslide bodies and collapse bodies that have already existed on the banks, and their changes will often further affect the stability of the reservoir bank slopes. Under the influence of external factors such as the sudden change of water level, some landslides that are in a relatively unstable state may collapse [3] . If the collapsing landslide has a large volume and a high speed, it will inevitably provoke a surge.
Scholars have conducted extensive research on the issue of landslide generated impulse waves and established various analytical methods and calculation formulas [4] . However, due to the complex interaction between the landslide body and the reservoir water and air, the landslide process is not only affected by slip speed, volume of instability, depth of water, but also by reservoir topography, reservoir width, duration of landslide process and length of landslide body [5] . Many of these factors affecting the generation and propagation of impulsive waves cannot be clearly determined and but only estimates, and the boundary conditions and initial conditions during are also complex, resulting in difficulties in simulation, and basically staying at the level of empirical calculation for a long time. In addition, in the canyon area, waves are also affected by factors such as obstacles on both sides of the river valley, refraction from the back and forth, and mutual interference of wave groups. This makes the prediction of landslides and tsunamis very complex and has not yet been solved satisfactorily [6] .
The most important aspect to truly simulate the generation and propagation of impulsive waves is to capture the free surface of the fluid. The currently used numerical simulation methods are mainly based on shallow water equations and Navier-Stokes(N-S) equations. The simulation method based on shallow water equation [7] can describe the process of generation and propagation of surges. However, there are many empirical parameters, and some parameters are difficult to determine in practical problems. The simulation method based on Navier-Stokes equations [8] requires less empirical parameters and has a great advantage in describing strong nonlinear waves. The SPH method (Smoothed Particle Hydrodynamics) is a meshless method that simulates continuous particles with smooth particles. It can describe the free surface well, facilitate the fluid solid coupling simulation, and can be used to simulate the large surge caused by high speed landslides. In recent years, it has been paid more and more attentions. For example, Capone [9] and Shi [10] simplified the landslide body into a block-shaped rigid body to simulate the generation and propagation of landslide generated impulse waves. Ataie-Ashtiani and Shobeyri [11] and Wang [12] considered the landslide body deformation and simulated the landslide impulsive waves under twodimensional simple boundary conditions. However, these works either do not take into account the deformation of the landslide body or the boundary conditions are relatively simple, so it is difficult to portray the characteristics of landslide impulsive waves with complex conditions, which limits its application in practical projects.
Considering the unique advantages of smooth particle dynamics method and discrete element method when dealing with the problem of large deformation of continuum and discrete media, this paper proposes a unified particle method that uses the coupling of them to deal with the interaction between discrete bodies and fluids at the macroscopic scale. This method was applied to simulate the two-dimensional wave generator experiment containing simple fluid-solid coupling to verify its effectiveness. Then the process of Huangping landslide generated impulse waves at the Xiluodu hydropower station of China was simulated and an effective way was provided to simulate the landslide and tsunamis.

Theoretical method
For fluid such as water, the N-S equation describing its motion is as follows: Where ρ is the density, v is the velocity vector, P is the pressure, g is the physical acceleration, and Θ is the viscous dissipation term. Reference Monaghan [13] gives the discrete form of the equation (1) from as follows: Where m is mass, subscripts a and b represent particles a and b, and ∇ a W ab is the first derivative of the kernel function. The viscosity term f ηa is used to characterize the viscous dissipative term Θ, using a simplified laminar viscosity solution proposed by Morris [14] . F a is the physical force received on a unit mass of fluid particle a.
The following commonly used weakly compressible assumption equations are used to determine the relationship between fluid pressure and density: Where ρ 0 is the reference density, and c 0 is the initial speed of sound. Generally for water, the constant γ=7.
The kernel function uses the cubic spline function which is widely used at present, and its expression is as follows: Where h is a smooth length, r is the distance between particle a and particle b, q = r/h = |r a -r b |/h, and α D is a regularization coefficient, which is equal to 10/(7πh 2 ) for 2D and 1/(πh 3 ) for 3D respectively.
For discrete solids such as rock and soil, the interparticle damping effect reflects the energy dissipation in the collision process. Considering the force of the particles under damping, the force can be expressed as [12] : Where the subscripts n and t indicate the normal and tangential directions, k is the elastic coefficient, λ and δ are the normal and tangential overlaps of the particles, c is the damping coefficient, and v sab is the relative slip speed of particle a and the particle b. n is the unit vector, and the normal and tangential elasticity coefficients are expressed by the following formula: Where E, G, ν, r are the elastic modulus, shear modulus, Poisson's ratio, and radius of the solid particles, respectively.
According to the principle of force synthesis, the final external force F a and the external moment T a of the solid particle a can be represented by the following formula: Where N is the number of particles b in the computing domain that are in contact with particle a, and the sliding friction between particles f s = min(f t , μf n ), where μ is the sliding friction coefficient.
According to Newton's second law of motion, the motion equation of solid particle a can be obtained as follows: v' a = F a /m a , w' a = T a /I a Where v' a and w' a are the acceleration and rotational acceleration of the particle a respectively. m a and I a are the mass and the moment of inertia.
When the fluid flows in a porous medium such as a discrete rock mass, energy is lost due to the viscous effect. According to Darcy's law, the expression of the viscous fluid in porous media is as follows: Δv = Q/A = KΔH/ηL (9) Where Δv is the relative flow velocity between the fluid and the porous medium; Q is the permeate flow; A is the flow area; K is the permeability coefficient, and its value should be determined by the properties of the solid medium, such as the material and pore size. The looser and more water-permeable solid medium K should have higher values, whereas the value of K is smaller; η is the dynamic viscosity coefficient of the fluid; L is the distance between the force sections. Substituting the head difference equation ΔH=ΔP/ρg into equation (9) results in: Δv = KΔP/ρgηL (10) Where ΔP is the pressure difference. Taking the minimum value of L in formula (10), the acceleration item of the fluid can be obtained as: ΔP/ρ = gηΔv/K (11) Multiplying the the acceleration by the density, and then integrating the volume of fluid particle group, the expression for the porous solid medium resistance f P of the fluid element can be obtained as follows: Where V is the volume of fluid particle group, the subscripts P, F represent solid media and fluid media.
Assuming that the fluid density and pressure are the same throughout the fluid, the expression of resistance f P and acceleration a F of solid medium to fluid element can be expressed as: f P = -ΔPV = -m P ΔP/ρ = -m P gηΔv/K (13) a F = -gηΔv/K (14) According to the principle of force and reaction force, the acceleration a P of fluid to solid medium can be expressed as: a P = -m F gηΔv/m P K (15) Using the particle approximation of the SPH to rewrite the acceleration a F and a P , the final equation for the interaction between the fluid particle and the discrete solid particle can be obtained as follows:

Validation
In this paper, the classical Scott Russell's wave generator experiment is used to verify the unified particle method [15] . In the simulation, the wave generated by the free fall of the slider at one end of the twodimensional tank will propagate to the other end of the tank. The simulation results were compared with the experimental results, and the validity of the method was verified by analyzing the characteristics of the flow field and waves. Ataie-Ashtiani and Shobeyri [11] proved that the length of the flume has little influence on the near-field characteristics in the experiment where the slider is vertically inserted into the water. Therefore, in order to improve the computational efficiency and obtain better visual effects in this paper, the length of the tank in the numerical model is shortened from 9.0 m to 2.0 m in the test tank. As shown in Fig.1, where d is depth of the sink.

Initial setup
The physical parameters and initial conditions of the numerical model are the same as those of the test model: the length of the tank is 2m and the width is 0.4 m. The three water depth conditions are 0.116 m, 0.21 m and 0.288 m respectively, and the water surface is 2 m long. The slider body length is 0.3 m, and the width is 0.4 m. The water density is 1000 kg/m 3 , and the slider body density is 816.24 kg/m 3 . The Poisson's ratio of slider is 0.3, the permeability coefficient K is 0.01, and the slider elastic modulus and shear modulus are 5.5×10 7 pa. The coefficient of collision recovery of the slider body and the sink is 0.25, and the friction coefficient of the sink is 0.25. The gravitational acceleration is in the +Y direction and the value is -9.81 m/s 2 , the initial velocity of the slider body into the water is 0, the distance between the particles is 0.005 m, and the time step is 1×10 -4 s.  [15] ; (b) simulated results. Table 1 gives a quantitative comparison of the experimental observations of the wave amplitude and the simulation results. The error between the wave amplitudes at different depths of water and the experimental results obtained by the simulation is almost 0, indicating that the method is effective.

Application
The unified particle method is properly validated with Scott Russell's wave generator experiment as shown before. In this section, the method is applied to study the real landslide and tsunamis. On the afternoon of July 27, 2013, landslide occurred in Huangping Village, Yongshan County, Yunnan province, China. The landslide site is shown in Fig.3. Sliding range is about 200 m in the direction of the river,110 m above the river level and 50 m above the Huangma road in some areas. The thickness of the sliding body is 3~5 m, and the landslide volume is about 120 thousand m 3 . The slippery body instantaneously entered the Jinsha River and formed waves more than 20 meters high. As a result, 12 people were missing on the other side (in the vicinity of the Yanjiao Township Pier under construction), 3 people were injured, and 2 speedboats were damaged [16,17] .

Numerical model and parameters
The sediment in the river bottom is reduced to water, and the unsteady shoreline rock and soil are regarded as landslide boundaries. Omitting the influence of air, the 2D numerical calculation model for establishing the Huangping landslide area is shown in Fig.4. Among them, the model consists of three parts: the sliding body, the water body, and the boundary. In order to observe the details of the landslide generated impulse waves, the monitoring area A where the landslide enters the water and the monitoring area B of the opposite bank are shown in dotted lines in the fig.4. Here, the surface width of the river is 745 m and the depth of the water is 130 m. The unstable landslide body is shown in Fig.4. The density of the landslide body particles is 2100 kg/m 3 , and the elastic modulus and shear modulus are 5.5×10 7 pa. The Poisson's ratio is 0.3, the permeability coefficient K is 0.0001, the coefficient of collision recovery between the sliding body particles and the water tank and itself is 0.25, the friction coefficient between the sliding body particle and the landslide boundary and itself is 0.38, the acceleration of gravity is -9.81 m/s 2 in the +Z direction, the particle spacing is 1 m, and the time step is 1×10 -6 s.

Results and discussions
The whole process of the sliding body swell includes five typical states: the falling of the landslide, the swell of the sliding body into the water, the propagation of the wave to the opposite bank, and the rising and falling of the wave in the opposite bank. The key moments are shown in Fig.5, where the red color and large arrows mean that the speed is greater. During the entire dynamic process of the landslide entering the water, the velocity of the water near the falling point and the surge front is relatively large. When the Huangping landslide occurred, the water level was 554 m, and the lowest unstable slippery body was in contact with water. Therefore, there was no obvious splash phenomenon when the slippery body impacted the water body, as shown in Fig.5(a)i. As the sliding body continuously slides into the water, the water level near the water fall rapidly increases until the sliding body almost completely enters the water. The water body is discharged to the left and generates a maximum surge of 24 m and propagates to the opposite shore, as shown in Fig.5(a)ii. As the swell spreads to the opposite shore, the body of water in front of the swell wave begins to climb, as shown in Fig.5(b)i-iii. After the surge reaches the maximum point of climb, the water body will fall back, and the potential energy will be converted into kinetic energy during the fall process, which will increase the speed of the water body, as shown in Fig.5 The only monitored impulse wave height at the local is an estimation value of 20 m. So the simulation value is larger than the estimated value of about 16.7%. However, as far as experience is concerned, the estimated 20 m is generally about 25 m. Therefore, the error between the numerical simulation results and the actual maximum wave height should be much smaller than 16.7%. This means that the accuracy and credibility of the method are relatively high.

Conclusions
In the current work, a unified particle method is presented to simulate the fluid-solid coupling problem in broad range of scales. In the first part of this paper, the Scott Russell wave generator experiment is simulated by using the unified particle method and compared with experimental results. It was found that the simulation results were in good agreement with the experimental results. The wave amplitude under different water depth conditions is within the experimental error range, that means the error is equal to 0, which indicates that the method is feasible and effective. In the second part of this paper, the 2D process of Huangping landslide generated impulse waves was further simulated, and the maximum wave height was compared with local observations. The results show that the numerical simulation results are close to the local observations.
We can clearly see the impact of the landslide body on water, the generation and propagation of waves, and the deformation process of the bulk bodies and water. Therefore, the unified particle method can provide an effective means for accurately and conveniently obtaining details in the process of landslide and tsunamis.