Numerical Runout Modeling Analysis of the Loess Landslide at Yining, Xinjiang, China

The Panjinbulake loess landslide is located in the western part of the Loess Plateau, in Yining County, Xinjiang, China. It is characterized by its long runout and rapid speed. Based on a field geological survey and laboratory test data, we used the DAN-W dynamic numerical simulation software (Dynamic Analysis Of Landslides, Release 10, O. Hungr Geotechnical Research Inc, West Vancouver, BC, Canada) and multiple sets of rheological models to simulate the whole process of landslide movement. The best rheological groups of the features of the loess landslide process were obtained by applying the Voellmy rheological model in the debris flow area and applying the Frictional rheological model in the sliding source area and accumulation area. We calculated motion features indicating that the landslide movement duration was 22 s, the maximum movement speed was 20.5 m/s, and the average thickness of the accumulation body reached 5.5 m. The total accumulation volume, the initial slide volume and the long runout distance were consistent with the actual situation. In addition, the potential secondary disaster was evaluated. The results show that the DAN-W software and related model parameters can accurately simulate and predict the dynamic hazardous effects of high-speed and long runout landslides. Together, these predictions could help local authorities make the best hazard reduction measures and to promote local development.


Introduction
Landslides are among the most destructive geological disasters with features of rapid speed, long runout distance, and entrainment effect [1][2][3]. Catastrophic landslide events are often triggered by heavy rainfall, earthquake, and engineering activities [4][5][6]. According to the spatial characteristics and trajectory of the sliding body, the entire process of landslide movement is mainly divided into three stages: The starting stage at the slide source area, the propagating stage, and the deposition stage [7][8][9][10]. In addition, entrainment, base liquefaction, and air cushioning occur during the landslide movement [11][12][13]. To reduce the landslide hazard loss, risk assessment is often requested [14][15][16]. The landslide runout analysis is a very effective method to assess a landslide hazard [17]. Landslide runout analysis involves two aspects: The simulation of previous landslides and the prediction of potential landslides [18]. Runout analysis could be used to design remedial engineering measures, such as barricades and berms [17]. The maximum runout distance, propagation velocity, and the deposit thickness, and provision of the basis for the design of remedial engineering measures, are obtained by landslide runout analysis.

Figure 1.
Location map of the Panjinbulake loess landslide (a is modified from [33]; b is modified from [48]). (a) The distribution of loess in China. (b) The distribution of loess in the study area.

Site Overview
The Panjinbulake loess landslide is located in the Karayagaq Township, Yining County, Xinjiang. The coordinates of the central point of the landslide are 81°30′32″ E and 44°11′48″ N, 53 km from downtown Yining City, and 43 km from Sailimu Lake (Figure 1b).
Among the landslide geological hazards in the Piliqinghe basin, most are high speed and long runout loess landslides, which pose a huge hazard to local agricultural and livestock production. Among these landslides, the Kezileisai landslide buried 402 cattle and 5 sheep, and the direct economic losses reached 490,600 Yuan (¥). There were some similar landslides in the Piliqinghe Basin ( Figure 2b). The Panjinbulake loess landslide is a typical high speed and long runout loess landslide that occurred recently in the region. Through analysis of the induced factors and dynamic effects of the landslide, this study provides a reference to the landslide dynamic hazards in the loess area. Figure 1. Location map of the Panjinbulake loess landslide (a is modified from [33]; b is modified from [48]). (a) The distribution of loess in China. (b) The distribution of loess in the study area.
In the study, we studied the characteristics of the Panjinbulake loess landslide through a field geological survey and aerial image analysis using drones. We used the landslide dynamics model DAN-W and multiple sets of rheological models to calculate the dynamic characteristics of this landslide. By using the different rheological models to simulate the different stages of the loess landslide (i.e., triggering in the sliding source area, propagation in the debris flow channel area, and deposition in the accumulation area), the best rheological model groups and parameters were obtained to improve the accuracy to analyze the loess dynamic characteristic. The potential secondary failure of the landslide was evaluated. This study could offer a basis to predict the potential landslide runout distance and define the hazard area, make necessary measures to prevent landslide induced damages (e.g., engineering measures, landslide early warning systems, and emergency response), and to favour local development.

Site Overview
The Panjinbulake loess landslide is located in the Karayagaq Township, Yining County, Xinjiang. The coordinates of the central point of the landslide are 81 • 30 32" E and 44 • 11 48" N, 53 km from downtown Yining City, and 43 km from Sailimu Lake (Figure 1b).
Among the landslide geological hazards in the Piliqinghe basin, most are high speed and long runout loess landslides, which pose a huge hazard to local agricultural and livestock production. Among these landslides, the Kezileisai landslide buried 402 cattle and 5 sheep, and the direct economic losses reached 490,600 Yuan (¥). There were some similar landslides in the Piliqinghe Basin (Figure 2b). The Panjinbulake loess landslide is a typical high speed and long runout loess landslide that occurred recently in the region. Through analysis of the induced factors and dynamic effects of the landslide, this study provides a reference to the landslide dynamic hazards in the loess area.

Post-Failure Behavior
The triggering of landslides was mainly from snow infiltration, which turn into long runout and rapid landslides, constitutes a typical disaster model in the loess area. The Panjinbulake landslide belongs to this typical geologic hazard model. Following the instability failure of the landslide, the front loess main body slipped due to the river flushing action. Pore-water pressure increased and soil saturation during the sliding process because of the snow infiltration. The landslide was transformed into debris flow, showing a flow state. Then, the landslide volume increased gradually with entrainment effects and flowed into the Pliqinghe Gully. Finally, the landslide struck the opposite side of the mountain and stopped.

Landslide Influential Factors
According to the field investigation, the occurrence of the Panjinbulake landslide was caused by a combination of factors, including geological structure and formation lithology, topography and hydrogeological condition.

Geological Structure and Formation Lithology
The landslide research area belongs to the western part of the Yili Valley. The terrain is generally high in the north and low in the south. It is gradually inclined from the northeast to the southwest.

Post-Failure Behavior
The triggering of landslides was mainly from snow infiltration, which turn into long runout and rapid landslides, constitutes a typical disaster model in the loess area. The Panjinbulake landslide belongs to this typical geologic hazard model. Following the instability failure of the landslide, the front loess main body slipped due to the river flushing action. Pore-water pressure increased and soil saturation during the sliding process because of the snow infiltration. The landslide was transformed into debris flow, showing a flow state. Then, the landslide volume increased gradually with entrainment effects and flowed into the Pliqinghe Gully. Finally, the landslide struck the opposite side of the mountain and stopped.

Landslide Influential Factors
According to the field investigation, the occurrence of the Panjinbulake landslide was caused by a combination of factors, including geological structure and formation lithology, topography and hydrogeological condition. The landslide research area belongs to the western part of the Yili Valley. The terrain is generally high in the north and low in the south. It is gradually inclined from the northeast to the southwest. The elevation of the area is between 620 and 3700 m. It is a block-like eroded and uplifted mountain, covered with gravel and loess layers, showing a low mountain grassland landscape. The landslide area is located on the southwestern side of the West Tianshan Youdi trough fold belt in the southwestern Tianshan fold of the Tianshan-Xing'an trough fold area. It belongs to the junction of the Boroconu Mountain Complex Anticline and Yili Block and is located 2.7 km south of the Nalati deep fault zone (Figure 2a). The rocks in the study area are mainly Ishikirik group tuff, tuff lava, and gray-green coarse sandstone in the Carboniferous system and basalt in the Dahala Junshan Formation of the Carboniferous system. The surface layer is the Quaternary Holocene loess, with well-developed joints (Figure 2b). The structure was relatively loose, and the wormhole, large void structure, belonging to low-plastic silt. The silt (0.075-0.005 mm) content of the loess in the Piliqinghe area is high, reaching 69.8-86.0%; the fine sand (0.25-0.075 mm) content is 3.7-18.0%; and the clay (<0.005 mm) content is 10.3-12.2%. As a result, the loess expands and collapses after encountering water, and it is prone to motion liquefaction under certain static or dynamic water conditions, which provide good source conditions for landslides.

Topography
The Panjinbulake loess landslide was located on the south bank of the river and had long been subjected to the lateral erosion of the river, resulting in good conditions for the landslide front to be in the air. The hillslope was steep, with a slope of 40 • (Figure 3a). Corresponding tensile stress condition occurred near the top of the mountain, and the cracks at the trailing edge of the landslide were gradually enlarged (Figure 3b). According to the Google Earth remote sensing image map from 18 May 2013, the front edge of the landslide had slipped. The sliding volume was about 9000 m 3 (Figure 3a), and several tensile cracks appeared on the trailing edge. Consequently, a steep ridge of up to 1 m was formed (Figure 3c,d), which provided good topographic conditions for loess deformation and stress relief.

Hydrogeological Condition
The landslide area belongs to the temperate continental semi-arid climate. The average annual precipitation (for the period of 2011-2016) is 330.6 mm. Precipitation is highest from March to July, during which the monthly rainfall exceeds 30 mm, accounting for 52.5% of the annual rainfall ( Figure 2c). Snowfall mainly occurs from October to the following March. The snowfall thickness can reach 94 mm per month. The fissure water inside the slope is frozen, and the vertical joints and cracks become enlarged due to the frost heaving action. From mid-March, the temperature rises, the snow that covered the surface begins to melt, infiltrating the cracks and joints and forming a certain transient water pressure and transient saturation zone in the surface layer of the slope. This results in a decrease in the anti-sliding force of the slope, thus inducing landslides. The rising water level of the river also causes the hydraulic gradient inside the slope to drop. Ice and snow meltwater can also be stored in the mountain for a long time and can continue to increase the slope sliding force and accelerate the formation of the potential slip surface. In addition, fissures are relatively developed and accumulated in the bedrock, which is exposed in the landslide. There is a large amount of ice and snow meltwater, which readily forms a "pipeline" channel that is in contact with the surface of the Quaternary aeolian loess and is discharged outward in the form of a spring. The flow volume of a spring was measured to be 5 × 10 −5 m 3 /s. The mineralization of water is less than 1.0 g/L. Based on the soil test, the saturation of the soil is 85.7-91.2%, which shows that it is a very wet sliding body. The natural moisture content of the soil was 17.2-20.4%, the plastic limit was 15.8%, and the liquid limit was 27.6%. According to the measured data of groundwater level, the groundwater in the hill is shallow and buried in the range

Hydrogeological Condition
The landslide area belongs to the temperate continental semi-arid climate. The average annual precipitation (for the period of 2011-2016) is 330.6 mm. Precipitation is highest from March to July, during which the monthly rainfall exceeds 30 mm, accounting for 52.5% of the annual rainfall ( Figure  2c). Snowfall mainly occurs from October to the following March. The snowfall thickness can reach 94 mm per month. The fissure water inside the slope is frozen, and the vertical joints and cracks become enlarged due to the frost heaving action. From mid-March, the temperature rises, the snow that covered the surface begins to melt, infiltrating the cracks and joints and forming a certain transient water pressure and transient saturation zone in the surface layer of the slope. This results in a decrease in the anti-sliding force of the slope, thus inducing landslides. The rising water level of the river also causes the hydraulic gradient inside the slope to drop. Ice and snow meltwater can also be stored in the mountain for a long time and can continue to increase the slope sliding force and accelerate the formation of the potential slip surface. In addition, fissures are relatively developed and accumulated in the bedrock, which is exposed in the landslide. There is a large amount of ice and snow meltwater, which readily forms a "pipeline" channel that is in contact with the surface of the Quaternary aeolian loess and is discharged outward in the form of a spring. The flow volume of a spring was measured to be 5 × 10 −5 m 3 /s. The mineralization of water is less than 1.0 g/L. Based on the soil test, the saturation of the soil is 85.7-91.2%, which shows that it is a very wet sliding body. The natural moisture content of the soil was 17.2-20.4%, the plastic limit was 15.8%, and the liquid limit was 27.6%. According to the measured data of groundwater level, the groundwater in the hill is shallow and buried in the range of 0-15 m. These provide good hydrological conditions for inducing Qualitative analysis of the disaster-formation mechanism of the Panjinbulake landslide was conducted based on a field geological survey and remote sensing satellite images, but this was far from adequate for geological disaster prevention and control. Instead, quantitative analysis methods are required to fully investigate landslide movement and predict the secondary disaster, which can be explored using the landslide dynamic analysis software DAN-W (see Section 5).

Basic Characteristics of Landslides and Hazard Zoning
The landslide had a long-strip shape ( Figure 4a). The slope before the landslide was close to 40 • , and the main slip direction was N 69 • E. The elevation of the trailing edge of the landslide was about 1280 m, the elevation of the landslide shearing edge was about 1190 m, and the horizontal distance reached 375 m (Figure 4b). Due to rainfall and snow infiltration, pore water pressure and soil saturation increased during the sliding process. The landslide was transformed into debris flow, showing a flow state. Based on the information obtained from the unmanned aerial vehicle data and remote sensing, geological field surveys, the landslide can be divided into the sliding source area, debris flow area, and accumulation area ( Figure 5).

Theoretical Basis
DAN-W is numerical simulation software developed by Hungr to simulate the whole process of landslide movement and to study the dynamics of landslides [23]. The 3D numerical model was set up according to the two-dimensional simulation conditions provided by the calculation profile in Figure 3b. Based on the aerial views, the path widths of landslide were confirmed. In DAN-W, the

Theoretical Basis
DAN-W is numerical simulation software developed by Hungr to simulate the whole process of landslide movement and to study the dynamics of landslides [23]. The 3D numerical model was set up according to the two-dimensional simulation conditions provided by the calculation profile in Figure 3b. Based on the aerial views, the path widths of landslide were confirmed. In DAN-W, the

Sliding Source Area (Area I)
The sliding body was mainly Quaternary loess. There were many cracks in the upper part of the trailing edge, the crack width ranged from 23 to 54 cm, and there was an unstable body in the upper part of the trailing edge. Its volume reached 1.2 × 10 4 m 3 . The average width of the sliding source area was 73 m. The area was about 1.2 × 10 4 m 2 , the average thickness of the sliding body was 4-7 m, and the volume was about 5.0 × 10 4 m 3 .

Debris Flow Area (Area II)
The debris flow area presented a long and narrow shape, and bedrock appeared in the upper part of the debris flow area. The bedrock surface had signs of scratch, and there was an exposed water head. There was also a water head on the western side of the sliding source area. This might suggest that there was sufficient groundwater in the area that promoted the sliding of the landslide. In this area, the volume of the sliding body increased due to the entrainment. The area was 6.8 × 10 3 m 2 , the average thickness of the sliding body was 3-6 m, and the volume was about 25.0 × 10 4 m 3 .

Accumulation Area (Area III)
The accumulation area had a fan shape, with a length of up to 110 m along the sliding direction and a maximum width of 100 m in the vertical sliding direction. The lithology of the accumulation area is dominated by Quaternary loess, which contains moderately weathered coarse sandstone scraped off the opposite mountain. The landslide struck the opposite side of the mountain and accumulated in the Piliqinghe Gully. The area was 1.46 × 10 4 m 2 , the thickness of the sliding body was 3-6 m, and the volume was about 6 × 10 4 m 3 .

Theoretical Basis
DAN-W is numerical simulation software developed by Hungr to simulate the whole process of landslide movement and to study the dynamics of landslides [23]. The 3D numerical model was set up according to the two-dimensional simulation conditions provided by the calculation profile in Figure 3b. Based on the aerial views, the path widths of landslide were confirmed. In DAN-W, the Lagrangian analytical solution of the Saint-Venant equation is mainly used to treat the sliding body with the rheological features that are formed by a combination of several blocks with certain materials ( Figure 6). In the curve coordinates, the corresponding physical equations and equilibrium equations are established for each block ( Figure 6), as in Equations (1)-(7) [23].
Here, F is the sliding force (N); γ is the unit weight (KN/m 3 ); H is the block height (m); B is the block width (m); ds is the nominal length of the block(m); α is the slope foot ( • ); P is the internal tangential pressure (N); and T is the base resistance (N); i is the block index.
Here, V is the new speed when sliding body movement. The new velocity at the end of a time step is obtained from the old velocity, v (m/s); g is the gravitational acceleration(m/s 2 ); ∆t is the time step interval(s); M is momentum flux; and the other parameters are the same as in Equation (1).
Here, h is the average depth of the slip mass; j is block boundary index; i is block index; S is the curve displacement (m); and the other parameters are the same as in Equation (1).
Here, V is the entire volume of the loess landslide deposits(m 3 ); V R is the volume of the initial landslide (m 3 ); V point is the volume of the unstable body (m 3 ); Y is the yield rate; L is the length of the i block; and i is the block index.
Water 2019, 11, x; doi: FOR PEER REVIEW www.mdpi.com/journal/water coefficient (f) is a parameter that should be determined using the Voellmy model. The friction coefficient was modified by the pore pressure and could reach much smaller values when the path material shows wet features. Frictional model: Assume that the flow of the sliding body is controlled by the effective normal stress acting on each block. The expression of resistance τ is as follows where γ is the pore pressure coefficient (specifically, the ratio of pore pressure to total stress); φ is the internal friction angle; σ is the total stress perpendicular to the direction of the sliding path; and τ is the resistance at the bottom of the sliding body.  (1)-(7) (It is modified from [23]).

Model Selection
In the DAN-W software, the accuracy of the calculation result depends on three important factors: Sliding body motion trajectory, rheological model, and parameter selection. First, based on multi-period remote sensing images, aerial imagery of drones, and field geological surveys, topographic lines before and after landslides were determined ( Figure 6). Second, the Panjinbulake loess landslide was divided into the sliding source area, debris flow area, and accumulation area ( Figure 5). According to the hazard characteristics of different regions, it is critical to select suitable rheological models for different regions. Since the sliding source area started from the shearing exit, there are signs of scratch on the exposed area of the bedrock. According to the existing research results [23,[49][50][51], the Frictional model(F) was suitable for the sliding source area. The Frictional model and the Voellmy model(V) were used in the debris flow area and the accumulation area, respectively. According to the landslide path sequence, four sets of the rheological model combinations, Frictional-Frictional-Frictional, Frictional-Frictional-Voellmy, Frictional-Voellmy-Voellmy, and Frictional-Voellmy-Frictioanl, were used to simulate the dynamic hazard effects of the  (1)-(7) (It is modified from [23]).
Momentum and mass during the entrainment of the path material could influence landslide kinematics. To describe the entrainment process quantitatively, an entrainment ratio (ER) could be offered to calculate the increase of the landslide volume for a specific entrainment zone in the DAN model [49].
where V E (i.e., V Entrained ) is the volume of the entrained path material (m 3 ); V Fragmented is the volume of the fragmented material in the sliding source area(m 3 ). V R is the volume of the initial loess landslide (m 3 ); and F F is the fractional amount of volume expansion due to fragmentation (0.25). The entire volume of the loess landslide deposits is equal to V R (1 + F F ) + V E [40]. In this study, V R equals 5.0 × 10 4 m 3 and V E equals 25.0 × 10 4 m 3 . The length of the entrainment area was approximately 240 m. To simulate the phenomenon of entrainment, an ER equal to 4.0 was used in the DAN model of the loess landslide. According to the pore-water pressure increased and soil saturation during the sliding process because of the snow infiltration. The landslide was transformed into debris flow, showing a flow state, so the scraping volume was huge. The movement speed of the sliding body and the thickness of the landslide accumulation body are calculated using Equations (1)- (5). In addition, the amount of resistance encountered during the movement of the sliding body is determined by different types of rheological models. In the DAN-W software, the resistance is mainly controlled by different base rheological models. DAN-W provides a range of rheological models. According to the existing research results and the trial-error method [50,51], the Voellmy model (V) and the Frictional model (F) are more suitable for landslide dynamic hazard research. The Frictional model is mainly used for landslides when the particle sizes of the residual body are large. The Frictional model is also used for mountains with open hillside cracks where the turbulent flow is not developed. The Voellmy model is suitable for the simulation of a landslide with fractured particles where there is a visible liquidized layer in the sliding mass. From Equation, it is evident that the rheological model is proportional to the velocity of the sliding body, so it could simulate the energy damage of the turbulent flow. This was caused by the liquefied material that has high moisture content, including the loose soil covering the flow path and a spring appearing in the path. This opinion has been accepted by Geotechnical Engineering Office (GEO) of Hong Kong [14,52].
Voellmy model: The expression of base resistance is as follows where f is the friction coefficient of the sliding body, σ is total stress perpendicular to the direction of the sliding path, γ is the material unit weight, ν is the moving speed of the sliding body, ξ is the turbulence coefficient, and τ is the resistance at the bottom of the sliding body. The constant friction coefficient (f) is a parameter that should be determined using the Voellmy model. The friction coefficient was modified by the pore pressure and could reach much smaller values when the path material shows wet features.
Frictional model: Assume that the flow of the sliding body is controlled by the effective normal stress acting on each block. The expression of resistance τ is as follows where γ µ is the pore pressure coefficient (specifically, the ratio of pore pressure to total stress); ϕ is the internal friction angle; σ is the total stress perpendicular to the direction of the sliding path; and τ is the resistance at the bottom of the sliding body.

Model Selection
In the DAN-W software, the accuracy of the calculation result depends on three important factors: Sliding body motion trajectory, rheological model, and parameter selection. First, based on multi-period remote sensing images, aerial imagery of drones, and field geological surveys, topographic lines before and after landslides were determined ( Figure 6). Second, the Panjinbulake loess landslide was divided into the sliding source area, debris flow area, and accumulation area ( Figure 5). According to the hazard characteristics of different regions, it is critical to select suitable rheological models for different regions. Since the sliding source area started from the shearing exit, there are signs of scratch on the exposed area of the bedrock. According to the existing research results [23,[49][50][51], the Frictional model(F) was suitable for the sliding source area. The Frictional model and the Voellmy model(V) were used in the debris flow area and the accumulation area, respectively. According to the landslide path sequence, four sets of the rheological model combinations, Frictional-Frictional-Frictional, Frictional-Frictional-Voellmy, Frictional-Voellmy-Voellmy, and Frictional-Voellmy-Frictioanl, were used to simulate the dynamic hazard effects of the Panjinbulake loess landslide, so as to select the most suitable rheological model combination to simulate the movement process of the loess landslide (see Table 1). The simulation results calculated from a combination of the four rheological models were compared to the features of the actual loess landslide. According to the results, we could find that the Frictional-Voellmy-Frictional model fit the above characteristics. From Equation, the Voellmy model is proportional to the velocity of the sliding body, and it could simulate the huge scraping force which removed the surface loess soil. Thus, the Voellmy rheology model had a better fit than the Frictional model in the debris flow area. Finally, we intended to simulate the movement characteristics of the Panjinbulake loess landslide using the Frictional-Voellmy-Frictional rheological model.

Parameter Selection
The dynamics of high speed and long runout landslides have been studied by researchers. The analysis of landslide dynamics depends, to a large extent, on the choice of parameters and the knowledge level of the author [23]. In this paper, the simulation parameters of the Panjinbulake loess landslide were mainly obtained by the field survey data and the existing research results [23,[49][50][51]. For the Voellmy model, the main parameters were ξ = 400 m/s 2 (the software provides a range of 200-500 m/s 2 ) and f = 0.05. When rain and snow melt water infiltrates, the groundwater level rises, and pore water pressure rises. The sliding body is close to the flow state, so the friction coefficient decreases gradually. The influence of underground groundwater on sliding body motion is realized by changing the friction coefficient in DAN-W (see Section 7.1). For the frictional model, the dynamic friction angle ϕ b was set to 19 • according the literature [23,[49][50][51]. Due to the infiltration of ice and snow meltwater, the excess pore water pressure increases, and Ru was set to 0.7. Finally, according to the indoor geotechnical test and the engineering analogy method, the typical strength testing index ϕi was set to 20 • and the unit weight (γ) was set to 18 KN/m 3 . As shown in Table 2, based on the trial and error method and the existing research results, these rheological model combinations and parameters were used to simulate the dynamic hazard effects of the Panjinbulake loess landslide.

Speed Analysis
Using the DAN-W software and the Frictional-Voellmy-Frictional model, the total time of the Panjinbulake loess landslide movement was 22 s. It is assumed that the speed of the landslide was 0 m/s when starting in the sliding source area. The sliding body started from the shear exit (i.e., the leading edge of landslide). Due to the steep slope of the landslide, the speed was very fast, reaching 5 m/s very quickly (Figure 7). Under the action of gravitational potential energy, the speed of the sliding body increased rapidly and reached 21.5 m/s at X = 200 m. This process lasted 7 s, and the average acceleration reached 2.3 m/s 2 . The sliding body moved from X = 200 m to X = 320 m in debris flow area. This area has exposed spring water ( Figure 5), which provides suitable hydrogeological conditions for high speed movement of the sliding body. Thus, the sliding body was always in a high speed state. This phase lasted 8 s. The average speed reached 20 m/s and the sliding body entered the river terrace. Due to the friction of sandy gravel and gradual decrease of the slope, the moving speed of the sliding body dropped sharply, and the speed relative to the opposite mountain was 17.5 m/s. Due to the blocking of the opposite mountain, the sliding body stopped at a horizontal distance of X = 370 m, and the average acceleration in the deceleration phase reached −3.6 m/s 2 , which lasted for 7 s.

Thickness Analysis of the Accumulation Body
As shown in Figure 8 The sliding body moved from X = 200 m to X = 320 m in debris flow area. This area has exposed spring water ( Figure 5), which provides suitable hydrogeological conditions for high speed movement of the sliding body. Thus, the sliding body was always in a high speed state. This phase lasted 8 s. The average speed reached 20 m/s and the sliding body entered the river terrace. Due to the friction of sandy gravel and gradual decrease of the slope, the moving speed of the sliding body dropped sharply, and the speed relative to the opposite mountain was 17.5 m/s. Due to the blocking of the opposite mountain, the sliding body stopped at a horizontal distance of X = 370 m, and the average acceleration in the deceleration phase reached −3.6 m/s 2 , which lasted for 7 s.

Thickness Analysis of the Accumulation Body
As shown in Figure 8, when the landslide sheared and started at X = 140 m, the initial volume of the landslide was 5.

Thickness Analysis of the Accumulation Body
As shown in Figure 8, when the landslide sheared and started at X = 140 m, the initial volume of the landslide was 5.62 × 10 4

Typical Point Velocity Analysis and Accumulation Body Thickness Analysis
Combined with the actual situation of the field geological survey, the four points with horizontal distances X = 140, 200, 250, and 350 m were selected as the typical points for analysis. The calculation results based on the DAN-W software and F-V-F rheological model are as follows.
(1) The point X = 140 m is located at the landslide shear outlet (Figure 4b). Within 0-10 s, the speed at this point increased linearly from 0 to 17.5 m/s (Figure 9a), which indicates that the acceleration of the sliding body during the starting process increased, reaching 1.75 m/s 2 . Also, typical negative terrain was present, and the lower part of the raised bedrock was exposed. This made it so that, after 12 s, the residual sliding body moved to the point where it was hindered by negative terrain with trough shape, the speed gradually attenuated, and the final thickness of the accumulation body at this point was 0.5 m (Figure 9b).
(2) The point X = 200 m is located in the debris flow area (Figure 4b). After the landslide slid for 7 s, it reached this point. Within 7-17.5 s, the speed at this location was always relatively faster, and the average speed reached 17.5 m/s (Figure 9c). In this area, the sliding body had a rapid speed and exposure to spring water, which provides conditions for the occurrence of a high speed and long runout landslide. At 8 s, the thickness of the accumulation body reached 4.53 m. After 17.5 s, the speed at this point gradually decreased, and the final thickness of the accumulation body reached 0.5 m (Figure 9d).
(3) The point X = 250 m is located at the foot of the slope (Figure 4b). After the sliding body slid for 9 s, it reached the foot of the slope. Due to the steep slope, the peak speed reached 19 m/s (Figure 9e), and the thickness of the accumulation body reached 5 m (Figure 9f). Scattering phenomenon was appeared at the foot of the slope where in the interval t = 20-24 s, and there was an increase of debris flow depth reaching a depth higher to the front passage. In addition, this was also a turning point for the speed of landslide movement. From this point on, the speed of the sliding body began to decrease due to the sudden slowing of the slope and the friction of the sliding body against the sandy gravel of the river terrace.
(4) The point X = 350 m is located at the foot of the opposite slope (Figure 4b). After the sliding body crossed the river terrace, it reached this point at 16 s; then, the speed dropped to 11 m/s (Figure 9g), and the thickness of the accumulation body reached 1.5 m (Figure 9h). From this point on, the sliding body began to hit the opposite side of the mountain. Because of this and even though a considerable amount of energy had been consumed, the sliding body continued to climb 20 m before stopping, indicating that a high-speed and long runout landslide has incredible energy and is severely devastating.

Sensitivity Analysis
The DAN-W dynamic model was preferably used to predict the dynamic characteristics of the landslide accumulation area in the Loess Plateau. However, based on a field survey and laboratory

Sensitivity Analysis
The DAN-W dynamic model was preferably used to predict the dynamic characteristics of the landslide accumulation area in the Loess Plateau. However, based on a field survey and laboratory tests, some parameters, e.g., the friction coefficient, were obtained by trial and error and the existing research results [23,[49][50][51]. There are many factors affecting the friction coefficient, including the rate of water content, terrain, and ground temperature, with the water content the primary factor. Therefore, it is difficult to provide the friction coefficient in a more efficient way. When using the DAN-W dynamic model to predict the risk assessment, the variability range of the parameters should be considered. In this paper, the influence of friction on the landslide is simulated using multiple sets of working conditions. According to the numerical values recommended by the kinematic model software and existing research results, the friction coefficient is divided into three groups. The first group has a low friction coefficient (i.e., 0.05, 0.1, and 0.2) which means high water content, the second group a moderate friction coefficient (i.e., 0.2, 0.3, and 0.4) which means moderate water content, and the third group a high friction coefficient (i.e., 0.5, 0.55, and 0.6) which means low water content. The calculation results are as follows ( Figure 10). group has a low friction coefficient (i.e.; 0.05, 0.1, and 0.2) which means high water content, the second group a moderate friction coefficient (i.e.; 0.2, 0.3, and 0.4) which means moderate water content, and the third group a high friction coefficient (i.e.; 0.5, 0.55, and 0.6) which means low water content. The calculation results are as follows (Figure 10). According to the numerical simulation results, it can be seen from Figure 10 that, if different friction coefficients are selected, the calculated moving distances are significantly different. For example, when a low friction coefficient (0.05-0.2) is selected, the moving distance of the sliding body is 350-375 m. When the moderate friction coefficient (0.2-0.4) is selected, the moving distance of the sliding body is 320-350 m. When the high friction coefficient (0.4-0.6) is selected, the moving distance of the sliding body is 280-300 m. According to the current geological survey in the field and the aerial view of the drone, the moving distance of the landslide reached a maximum of 366-375 m, which is in line with the calculation results based on the low friction coefficient. The low friction coefficient appears to be due to the saturation degree of its material, and there was sufficient groundwater in the slope area that promoted the sliding of the landslide.

Empirical-Statistical Model
The empirical formula, mainly based on simple geometric relations of landslides (Figure 11), is simple and effective in the prediction of long runout distance. Figure 11 reveals the geometric relationship between the landslide's apparent friction angle (i.e.; the angle between the trailing edge of the landslide and the farthest point of landslide movement), height difference, and motion distance [19]. Based on the concept of the apparent friction angle, Scheidegger proposed an empirical formula called the sled model to calculate the velocity of the sliding body [19]. The specific formula is as follows where is the velocity of the estimating point (m/s); g is the acceleration of gravity (m/s 2 ); and t is the rake ratio between the highest point of the rear edge and the estimating point of the sliding distance (dimensionless); H is the height difference from the highest point of the rear edge and the calculated velocity point of the sliding distance(m); L is the horizontal distance between the landslide trailing edge and the calculated velocity point (m). According to the numerical simulation results, it can be seen from Figure 10 that, if different friction coefficients are selected, the calculated moving distances are significantly different. For example, when a low friction coefficient (0.05-0.2) is selected, the moving distance of the sliding body is 350-375 m. When the moderate friction coefficient (0.2-0.4) is selected, the moving distance of the sliding body is 320-350 m. When the high friction coefficient (0.4-0.6) is selected, the moving distance of the sliding body is 280-300 m. According to the current geological survey in the field and the aerial view of the drone, the moving distance of the landslide reached a maximum of 366-375 m, which is in line with the calculation results based on the low friction coefficient. The low friction coefficient appears to be due to the saturation degree of its material, and there was sufficient groundwater in the slope area that promoted the sliding of the landslide.

Empirical-Statistical Model
The empirical formula, mainly based on simple geometric relations of landslides (Figure 11), is simple and effective in the prediction of long runout distance. Figure 11 reveals the geometric relationship between the landslide's apparent friction angle (i.e., the angle between the trailing edge of the landslide and the farthest point of landslide movement), height difference, and motion distance [19]. Based on the concept of the apparent friction angle, Scheidegger proposed an empirical formula called the sled model to calculate the velocity of the sliding body [19]. The specific formula is as follows where V is the velocity of the estimating point (m/s); g is the acceleration of gravity (m/s 2 ); and t is the rake ratio between the highest point of the rear edge and the estimating point of the sliding distance (dimensionless); H is the height difference from the highest point of the rear edge and the calculated velocity point of the sliding distance(m); L is the horizontal distance between the landslide trailing edge and the calculated velocity point (m).
Water 2018, 10, x FOR PEER REVIEW 16 of 21 Figure 11. Sketch of the Empirical-statistical model (It is modified from [19]).
According to the calculation of the sled model (Figure 12), the maximum speed of the sliding body was 28 m/s, which occurred near the horizontal distance of 250 m. Similar to the calculation results of the DAN model, the sled model showed that the slide velocity increased sharply after the landslide occurred, and the speed decreased significantly when the slide moved to the front of the road and eventually struck the opposite side of the mountain. The maximum speed obtained by the sled model is far greater than that of the DAN-W dynamic model. Instead of taking the dynamic characteristic of the landslide into account, such as erosion and entrainment, the sled model only gives a preliminary description of the process of landslide movement variation. It could be seen that the calculation results of the DAN model are more accurate. While the concept of the sled model and the apparent friction angle tends to be conservative for landslide hazard prediction, they still comprise a qualitative and effective way to predict disasters. The landslide under study is located in Piliqinghe Basin, located in the western part of the Loess Plateau and is part of the "Belt and Road" area. The location has many potential loess landslides, all of which pose a threat to agricultural production. It is especially important to research the prediction of potential landslide disasters, which is of great benefit to disaster prevention and devising mitigation measures. During the period 2017-2018, our team carried out field geological survey work Figure 11. Sketch of the Empirical-statistical model (It is modified from [19]).
According to the calculation of the sled model (Figure 12), the maximum speed of the sliding body was 28 m/s, which occurred near the horizontal distance of 250 m. Similar to the calculation results of the DAN model, the sled model showed that the slide velocity increased sharply after the landslide occurred, and the speed decreased significantly when the slide moved to the front of the road and eventually struck the opposite side of the mountain. The maximum speed obtained by the sled model is far greater than that of the DAN-W dynamic model. Instead of taking the dynamic characteristic of the landslide into account, such as erosion and entrainment, the sled model only gives a preliminary description of the process of landslide movement variation. It could be seen that the calculation results of the DAN model are more accurate. While the concept of the sled model and the apparent friction angle tends to be conservative for landslide hazard prediction, they still comprise a qualitative and effective way to predict disasters.
Water 2018, 10, x FOR PEER REVIEW 16 of 21 Figure 11. Sketch of the Empirical-statistical model (It is modified from [19]).
According to the calculation of the sled model (Figure 12), the maximum speed of the sliding body was 28 m/s, which occurred near the horizontal distance of 250 m. Similar to the calculation results of the DAN model, the sled model showed that the slide velocity increased sharply after the landslide occurred, and the speed decreased significantly when the slide moved to the front of the road and eventually struck the opposite side of the mountain. The maximum speed obtained by the sled model is far greater than that of the DAN-W dynamic model. Instead of taking the dynamic characteristic of the landslide into account, such as erosion and entrainment, the sled model only gives a preliminary description of the process of landslide movement variation. It could be seen that the calculation results of the DAN model are more accurate. While the concept of the sled model and the apparent friction angle tends to be conservative for landslide hazard prediction, they still comprise a qualitative and effective way to predict disasters. The landslide under study is located in Piliqinghe Basin, located in the western part of the Loess Plateau and is part of the "Belt and Road" area. The location has many potential loess landslides, all of which pose a threat to agricultural production. It is especially important to research the prediction of potential landslide disasters, which is of great benefit to disaster prevention and devising mitigation measures. During the period 2017-2018, our team carried out field geological survey work The landslide under study is located in Piliqinghe Basin, located in the western part of the Loess Plateau and is part of the "Belt and Road" area. The location has many potential loess landslides, all of which pose a threat to agricultural production. It is especially important to research the prediction of potential landslide disasters, which is of great benefit to disaster prevention and devising mitigation measures. During the period 2017-2018, our team carried out field geological survey work in the area and counted 12 loess landslides that occurred. At the same time, the team measured the basic parameters of the landslide, and calculated the apparent friction angle of each landslide (Table 3, Figure 11). The statistical results indicate that the apparent friction angle of the loess landslide in this area is approximately 25 • . Based on the concept of the apparent friction angle, the farthest running distance of the landslide (i.e., L max ) can thus be calculated by the formula which is as follows where H max is the height difference from the highest point of the rear edge and the farthest point of the sliding distance (m). The empirical formula and the DAN model can be used to predict and analyze the moving distance of potential landslides from qualitative and quantitative aspects separately.

Evaluation of Landslide Residual Risk
Based on the field survey and the aerial images from unmanned aerial vehicles, multiple tensile cracks appeared in the upper part of the trailing edge of the landslide and formed an independent unstable block. In the case of rainfall and ice and of snow melt infiltration, the unstable block would be extremely easy to slide, which poses a threat to agricultural production and road operation. To avoid secondary harm, in this paper, we use the DAN model to predict and analyze the movement trend of the unstable body. According to the results of the field geological survey, the unstable block area ranges from 3 to 6 m, with a volume of nearly 1.2 × 10 4 m 3 .
The selected models and parameters are the same as in the previous landslide dynamic hazard analysis. Figure 13 shows the running velocity of the unstable body, the thickness of the deposit, and the predicted disaster threat zone. Since the slope of the unstable body surface is 30 • , which is relatively gentle, most of the sliding body is deposited in the slide-source area after the unstable body slides. The average thickness of the sliding body in the sliding source area reaches 5 m, the longest distance reaches 300 m, and the maximum moving speed reaches 18.5 m/s. According to the sled model, the farthest distance of the unstable body motion is 362 m. Combined with the calculation results of the DAN-W model and sled model, it can be concluded that, if the unstable landslide body starts, it will be pose a threaten to road operation. Owing to its rapid movement speed, the landside mass also threatens the safe production of grazing herds of animals. Setting up a warning sign around the landslide to warn herders to locate grazing far from the area is recommended. The local should set engineering measures (such as garbion) around the road to ensure traffic safety.

Conclusions
Based on the geological survey in the field, multi-period historical remote sensing images and aerial images of the drone, combined with the geological conditions of the study area, we analyzed the inducing factors and runout process of the Panjinbulake loess landslide and predicted the secondary disaster. Furthermore, the DAN-W dynamic model and a set of combined basal rheological models (Frictional-Voellmy-Frictional models) can suitably simulate the dynamic hazard effects of the Panjinbulake loess landslide. We analyzed the influence of the landslide movement speed, typical point velocity, accumulation body thickness, and friction coefficient. The simulation results showed that the duration of the Panjinbulake loess landslide was 22 s, the maximum speed was 20.5 m/s, and the maximum thickness of the accumulation body was 5.5 m, which is in line with the actual situation based on the field investigation. The basal rheological model combination and parameters obtained through trial and error can be used to simulate and predict the long runout distance of loess landslides, and it is necessary in strengthening the early identification and prevention of loess landslide hazards using multi-precision observation technology and numerical techniques.

Conclusions
Based on the geological survey in the field, multi-period historical remote sensing images and aerial images of the drone, combined with the geological conditions of the study area, we analyzed the inducing factors and runout process of the Panjinbulake loess landslide and predicted the secondary disaster. Furthermore, the DAN-W dynamic model and a set of combined basal rheological models (Frictional-Voellmy-Frictional models) can suitably simulate the dynamic hazard effects of the Panjinbulake loess landslide. We analyzed the influence of the landslide movement speed, typical point velocity, accumulation body thickness, and friction coefficient. The simulation results showed that the duration of the Panjinbulake loess landslide was 22 s, the maximum speed was 20.5 m/s, and the maximum thickness of the accumulation body was 5.5 m, which is in line with the actual situation based on the field investigation. The basal rheological model combination and parameters obtained through trial and error can be used to simulate and predict the long runout distance of loess landslides, and it is necessary in strengthening the early identification and prevention of loess landslide hazards using multi-precision observation technology and numerical techniques.