OPEN ACCESS
Understanding landslide thrust is the key basis for mitigating landslides. Since the spatial distribution of thrust cannot be acquired through traditional methods, there are conditions of some blindness and waste in actual treatment projects. A method is proposed to obtain the spatial distribution of thrust through using a typical accumulation landslide in Southwest China as an example, which is continuous homogeneous. Firstly the landslide is divided into grids with m rows and n columns. Secondly the value of thrust on each grid node is calculated through the MATLAB compiler based on the transfer coefficient method. Then the results are shown in the form of a threedimensional curved surface, and the distribution of thrust on any profile can be determined. The results indicate that landslide thrust on longitudinal profiles increases with the increase of the design safety factor, but the increasing range decreases. Furthermore, landslide thrust decreases with the decrease of angle of sliding surface or thickness. On cross profiles the distribution of thrust approximates broken lines. The spatial distribution visualization method of landslide thrust could provide a useful reference for practical engineering.
Accumulation landslide, Spatial distribution visualization, Thrust, Transfer coefficient method.
Landslides are a harmful and common physicogeological phenomenon all over the world, which not only pose a serious threat to the safety of life and property, destroying the natural scenery and landscapes, but also seriously affect the coordinated and sustainable socioeconomic development [13]. The stability analysis of landslides is related to the safety of human lives directly and environmental resources, so effective investigation and landslide treatment have attracted wide attention in recent decades [45].
Understanding landslide thrust is the key basis for designing landslide support [67]. In traditional investigations, a number of profiles along the main sliding direction and some vertical profiles were surveyed in detail. In general, calculating landslide thrust on these important profiles could meet the requirement of landslide control. But due to the varying thrust in different positions, designing landslide support according to thrust alone on these profiles could lead to waste and hazardous conditions [8]. It is necessary to study the spatial distribution of landslide thrust under limited conditions to conclude the weak position of landslides [9].
This article proposes a method to study the spatial distribution of thrust. Firstly topographic information of landslide surface and sliding surface is obtained through geological investigation, and the curved surface is fitted through SURFER software. Secondly the research area is divided into grids with m rows and n columns. Thirdly the thrust of each grid node is calculated automatically based on transfer coefficient method through MATLAB compiler. Finally the calculated result is shown in the form of threedimensional curved surface, where spatial distribution of landslide thrust could be observed visually. And on this basis, thrust on any profile can be sliced out of the 3D curved surface to analyze its distribution regularities.
The transfer coefficient method, with the simple force system, is considered as standard for geohazard investigation, which is widely adopted in engineering practice for landslide investigation, treatment and mechanism analysis in China [1012]. In its original principle, the sliding mass above the sliding surface or the sliding zone is divided into a number of vertical slices artificially. The curved sliding surface on the profiles is approximated to polylines. Relatively, in this article, the slices are divided in quantity in a fixed size automatically, and the thrust of each grid node is calculated through MATLAB compiler.
2.1 Grid generation
After the main slip direction and the range is confirmed, the landslide zone is divided into grids with m rows and n columns. The Y axis is opposite to the main slip direction. On longitudinal profiles of landslides, a number of slices are formed through drawing a group of parallel lines downwards vertically, which intersect with the landslide surface and the sliding surface. The bottom lines of slices are straight lines made by connecting adjacent points. When the number of slices is large enough, the polylines approximate the actual sliding curved lines. The schematic diagram of grid generation is shown in Figure 1.
Figure 1. Grid generation of landslides
2.2 Formula derivation
The transfer coefficient method is based on the limit equilibrium theory, and assumes that the moving direction from slice i–1 (higher) to slice i (lower) is parallel to the slip zone of slice i–1. Figure 2 shows the forces acting on the sliding slice, where W_{i} is the weight of the slice, Q_{i} is the horizontal component of external forces (potentially including earthquake force, anchor cable force, and surface load, etc.), V_{i}is the vertical component of external forces, U_{bi} is hydraulic pressure of the ground water, R_{i} is the resisting force acting on slice i, P_{i1} is thrust of slice i1, which is essentially interslice force between slice i1 and slice i, P_{i} is thrust of slice i, α_{i} is the angle of the sliding surface, and S_{i} is area of the sliding surface, with the grid spacing as the width.
Figure 2. Forces acting on the slice
Basing on the conditions of static equilibrium in a normal direction and a tangential direction of sliding surface in one slice, the following formulas can be derived.
$\begin{aligned} P_{i} &=\left(W_{i}+V_{i}\right) \sin \alpha_{i}+Q_{i} \cos \alpha_{i} \\ &\left\{c_{i} S_{i}+\tan \varphi_{i}\left[\left(W_{i}+V_{i}\right) \cos \alpha_{i}+Q_{i} \sin \alpha_{i}U_{b i}\right]\right\} / F_{s} \\ &+P_{i1}\left[\cos \left(\alpha_{i1}\alpha_{i}\right)\tan \varphi_{i} \sin \left(\alpha_{i1}\alpha_{i}\right) / F_{s}\right] \end{aligned}$(1)
That can be expressed as,
$P_{i}=T_{i}R_{i} / F_{s}+\psi_{i} P_{i1}$(2)
With
$\psi_{i}=\cos \left(\alpha_{i1}\alpha_{i}\right)\sin \left(\alpha_{i1}\alpha_{i}\right) \tan \varphi_{i} / F_{s}$(3)
$T_{i}=\left(W_{i}+V_{i}\right) \sin \alpha_{i}+Q_{i} \cos \alpha_{i}$(4)
$R_{i}=\left[\left(W_{i}+V_{i}\right) \cos \alpha_{i}U_{b i}Q_{i} \sin \alpha_{i}\right] \tan \varphi_{i}+c_{i} S_{i}$(5)
where T_{i} is the sliding force except thrust, ψ_{i} is the transfer coefficient.
Eq(1) leads to the nth thrust as follows
$P_{\mathrm{n}}=\left(\sum_{i=1}^{n1}\left(T_{i} \prod_{j=i}^{n1} \psi_{j}\right)+T_{n}\right)\left(\sum_{i=1}^{n1}\left(R_{i} \prod_{j=i}^{n1} \psi_{j}\right)+R_{n}\right) / F s$(6)
When thrust transfers from the highest to lowest, the last thrust should become 0 in the limit equilibrium state. So in Eq(6), when the nth thrust P_{n} equals to 0, the Fs is the safety factor at that time, which is derived as follows.
$F s=\frac{\sum_{i=1}^{n1}\left(R_{i} \cdot \prod_{j=i+1}^{n} \psi_{j}\right)+R_{n}}{\sum_{i=1}^{n1}\left(T_{i} \cdot \prod_{j=i+1}^{n} \psi_{j}\right)+T_{n}}$(7)
In practical engineering, the nth thrust is calculated according to the design safety factor, which is generally greater than the calculated safety factor, thus the calculated nth thrust is greater than 0.
A typical landslide in Southwestern China can be used as an example to illuminate the method of spatial distribution visualization landslide thrust.
3.1 Engineering geological conditions
The Jiuzi landslide traverses the borders of Yunnan and Guizhou provinces in Southwest China. The accumulation landslide mass lies in the landforms of highland, in the front of which the first terrace of the Qingshui River is formed through the erosive activity. Its main slip direction is roughly 145 ° SE. The landslide body is roughly 23.4 m thick and 770 m long, with a narrow upper part (approximately 80–120 m) and a wide lower part (approximately 600–700 m). The landslide covers an area of 36.5 ha. Its volume is estimated to be 8.5×10^{6} m^{3}.
In the geological investigation, comprehensive measures were carried out to clarify the geological structure and physicsmechanics characteristic of landslides, including engineering geological mapping, digging trenches and adits, drilling, and geophysical prospecting [13]. The landslide’s body consists of dolomite of Huanglong Formation, Carboniferous. And limestone of Baizuo Formation developed under the sliding surface. The rock mass has become compressed fragmentized by the action of the sliding friction, and it was cemented again by calcium with groundwater. The shape and location of the sliding surface could be inferred through drilling, trenches, adits, and geophysical prospecting. The surface boundary of the landslide was identified according to engineering geological mapping, and the distribution characteristic of rock, soil and vegetation at the site. The rock mass of the landslide accumulation can be regarded as continuous homogeneous due to its being loose and fractured. The unit weight γ is 26.3 kN/m^{3}, cohesion c is 150 kPa, and the internal friction angle φ is 28 °. The topographic map of the landslide with contour interval of 10 m is presented in Figure 3.
Figure 3. Topographic map of landslide
(I, highland; I1, landslides; II, first terrace of Qingshui River; Q^{del}, landslides accumulation; C_{2hn}, dolomite of Huanglong Formation, Carboniferous; C_{1b}, limestone of Baizuo Formation; C_{1y}, limestone of Yanguan Formation; 11',22'.., investigation profiles.)
3.2 Surface and sliding surface of landslide’s fitting
The main slip direction and the range of the landslide have been assessed according to the geological investigation. The quantity of horizontal coordinates and elevation data of the landslide surface can be extracted from CAD topographic map. The curved surface of the landslide was fitted through SURFER software, and then the surface was divided into grids with m rows and n columns.
Threedimensional coordinates of 41 points on the sliding surface were obtained through drilling. And some threedimensional coordinates of boundary points were extracted from a CAD topographic map. Then the sliding surface was fitted through SURFER software according to these coordinates. After that the sliding surface was divided into grids with m rows and n columns as well, following the same method as was used for the surface. In the rectangular coordinates, the values of m and n are both 160, and the distance between adjacent nodes is 5 m. In summary, there are 25,600 nodes in the research area.
The overlay results of the surface and the sliding surface are shown in Figure 4. The upper grid is the surface, and the lower is the sliding surface. The contour lines of 50 m are added to make the surfaces clearer.
Figure 4. Grids of surface and sliding surface
3.3 Calculation of 3D distribution of thrust
Coordinates of any grid node and mechanical parameters had been obtained through the above analysis. The thrust of any grid node acting on slices and the safety factor of any longitudinal profile can be calculated through Eq(6) and Eq(7). The value of thrust can be adjusted to correspond to the design safety factor. The design thrust increases with the design safety factor. The whole computational process was programed through MATLAB compiler, and the fixed MATLAB compiler can be applied in similar landslides engineering.
When the design safety factor is 1.2, 3D curved surfaces of thrust are shown in Figure 5. These four sub graphs are shown in different angles of view, with the format “ azimuth angle, elevation angle”. (a)~(d) are listed as (45, 25), (45, 25), (90, 0), (0, 0).
Figure 5. Thrust distribution when the safety factor is 1.2
When the design safety factor is 1.2, thrust increases first then decreases, and after that increases constantly along the slip direction from the higher to the lower, with the largest thrust in the front end of the landslide. So if the supporting structure is set there, the cost and the risk would both be high.
There is an obvious concave shape near the front, where the lower thrust is located. If the first line of antislide piles and retaining wall are set in this position, the cost of landslide control would be lower. And the position of the second line of the landslide support is reselected after the thrust of the first line is limited at 0.
It can be seen from the analysis that the spatial distribution visualization method can provide reference for setting up landslide support while at the same time reducing the cost.
Thrust on any profile can be sliced out from the 3D curved surface of thrust, and safety factors of longitudinal profile can be calculated through the programed MATLAB compiler as well.
4.1 Distribution regularities of safety factors
The stability level in landslides is different at different positions. When landslides begin to slide, the parts with poor stability slide firstly, and then followed by the parts with better stability, resulting in the complete instability in the end. It is helpful to analyze the safety factors of longitudinal profiles as well. The result is presented in Figure 6.
Figure 6. Distribution of safety factors
From the distribution of safety factors, we can see that some safety factors on the left side of landslides are less than 1, but the safety factors on the right side are mostly greater than 1, with some even greater than 2. Corresponding to the shape characteristic of landslides, a scarp develops on the left trailing edge of the landslide. So it can be presumed that this position has poor stability, which requires more attention in the project of control. If the landslide starts to slip, the left part would slip first. The distribution of safety factors could provide a reference for the global stability of landslides.
4.2 Distribution regularities of thrust on longitudinal profiles
The longest longitudinal profile is sliced out from the 3D thrust curved surface to analyze the thrust distribution, of which the position is located in x=15 m in this example. Its projecting position is presented in Figure 7.
Figure 7. Position of the profiles
Figure 8 shows the natural shape of this profile. The safety factor in a natural state is 1.174. The design thrust varies with the design safety factor. The curves of thrust are shown in Figure 9 when the safety factor of the profile increases from 1.1 to 1.5.
Figure 8. The natural shape of longest longitudinal profile (x=15 m)
Figure 9. Curves of thrust under different safety factors
Figure 9 shows that when the safety factor is less than 1.174, the thrust of the landslide front is 0, which means landslide support is not necessary in this position under the smaller design safety factor. The thrust increases as the design safety factors increase, which is in line with the common opinion. Figure 9 also shows a<b<c & d<e<f<g, so it can be concluded that increasing range decreases as the design safety factors increase.
In addition, there are two turning points on the curves. One is nearly at y=450 m, and another is nearly at y=120 m. Through contrasting the thrust distribution with the shape of the profile, we can find that the angle of the sliding surface becomes smaller at y=450 m nearby. It can be concluded that thrust decreases as the angle of the sliding surface gets smaller. The thickness of landslides decreases at y=120 m nearby. So we can infer that thrust decreases with the decrease of the landslide’s thickness as well.
4.3 Distribution regularities of thrust on cross profiles
Six cross profiles with the same spacing are sliced out from the 3D curved surface of thrust, of which positions are at y=110 m, y=210 m, y=310 m, y=410 m, y=510 m, y=610 m. The projecting positions of these profiles are presented in Figure 7 as well, and the distribution of thrust is shown in Figure 10.
Figure 10. Distribution of thrust on cross profiles
Figure 10 shows that all curves of thrust are high in the middle and low in both ends. For studying the distribution regularities on cross sections, the six curves are fitted with two commonly used functions, quadratic function curves and broken lines. The fitting results are shown in Table 1, where RSS is the residual sum of squares, and R^{2} is the correlation coefficient.
Table 1. Fitting results
Position 
Broken lines 
Quadratic function 

RSS 
R^{2} 
RSS 
R^{2} 

y=110 m 
5.48×10^{10} 
0.985 
2.40×10^{11} 
0.742 
y=210 m 
8.78×10^{10} 
0.946 
9.67×10^{10} 
0.776 
y=310 m 
1.10×10^{11} 
0.915 
9.57×10^{10} 
0.739 
y=410 m 
7.05×10^{10} 
0.9505 
3.07×10^{10} 
0.953 
y=510 m 
2.05×10^{10} 
0.9745 
4.05×10^{10} 
0.826 
y=610 m 
5.70×10^{8} 
0.9905 
3.39×10^{9} 
0.845 
According to the fitting results, the fitting effect of broken lines is better than quadratic function curves, meaning that the distribution of thrust is nearly broken lines. In the progress of landslide treatment, if the conditions of investigation are limited, the broken lines formed by connecting the maximum thrust with the thrust of landslide edge can be regarded as the thrust distribution.
This paper proposes a spatial distribution visualization method of landslide thrust through a case study of the Jiuzi landslide in Southwest China. This kind of 3D curved surface of thrust can provide a visual display which can help engineers to choose the optimal position of landslide support and analyze the stability of the landslide.
The thrust of any position on the landslide and safety factors of any profile can be obtained through this spatial distribution visualization method.
On the longitudinal profiles, thrust increases as the design safety factors increase, but the increasing range decreases. Thrust decreases with the decrease of angle of the sliding surface and thickness of landslide.
On the cross profiles, distribution of thrust shows high in the middle and low in both ends, which approximates broken lines. When the conditions of investigation are limited, the broken lines, formed by connecting the maximum thrust with the thrust of the landslide’s sides, could be regarded as the distribution of thrust.
The authors would like to Jilin Provincial Investigation Design & Research Institute of Water Conservancy & Hydropower, that provided investigation information for this project.
Wi 
Weight of slice, kN 

Q_{i} 
Horizontal component of external forces (including earthquake force, anchor cable force, and surface load, etc.) , kN 

V_{i} 
Vertical component of external forces (including earthquake force, anchor cable force, and surface load, etc.), kN 

U_{bi} 
Hydraulic pressure of the ground water, kN 

R_{i} 
Resisting force acting on slice i, kN 

T_{i} 
Sliding force except thrust, kN 

P_{i} 
Thrust of slice i, kN 

Greek symbols 

α_{i} 
Angle of sliding surface, ° 

S_{i} 
Area of sliding surface, m^{2} 

ψ_{i} 
Transfer coefficient. 

γ 
Unit weight, kN/m^{3} 

c 
Cohesion, kPa 

φ 
Internal friction angle, ° 
[1] E. Hoek and J. BRAY, Rock Slope Engineering, revised second edition, Publication of Institution of Mining & Metallurgy, 1977.
[2] G. X. Wang, “Key technique in landslide control and its handling measures,” Chinese Journal of Rock Mechanics and Engineering, vol. 24, no. 21, 3818382, 2005.
[3] M. Shi, J. Chen, D. Sun and C. Chen, “Hazard assessment of debris flows based on the catastrophe progression method: A case study from the Wudongde dam site,” International Journal of Heat & Technology, vol. 33, no. 4, 217220, 2015. DOI: 10.18280/ijht.330429.
[4] K. Y. Choi and R. W. M. Cheung, “Landslide disaster prevention and mitigation through works in Hong Kong,” Journal of Rock Mechanics and Geotechnical Engineering, vol. 5, no. 5, 354365, 2013. DOI: 10.1016/j.jrmge.2013.07.007.
[5] L. Fu and Y. C. Jin, “Investigation of nondeformable and deformable landslides using meshfree method,” Ocean Engineering, vol. 109, 192–206, 2015. DOI: 10.1016/j.oceaneng.2015.08.051.
[6] C. H. Liu, C. X. Chen and X. T. Feng, “Discussion on design method of prestressed cable for soil slope,” Rock and Soil Mechanics, vol. 27, no. 8, 13491352, 2006. DOI: 10.16285/j.rsm.2006.08.024.
[7] Y. H. Xia and S. W. Bai, “Application of transfer coefficient method to scheme design of landslide treatment with slope cutting,” Chinese Journal of Rock Mechanics and Engineering, vol. 27, supp. 1, 32813285, 2008.
[8] J. Zhang, R. L. Hu and Z. Q. Li, “Distribution laws of thrust of talus slide on antisliding piles,” Chinese Journal of Geotechnical Engineering, vol. 34, no. 11, 20052010, 2012.
[9] E. Ausilio, E. Conte and G. Dente, “Stability analysis of slopes reinforced with piles,” Computers and Geotechnics, vol. 28, no. 8, 591611, 2001. DOI: 10.1016/S0266352X(01)000131.
[10] R. N. Bi, E. Dominik, W. Xiang, J. Joachim, S. Markus and J. W. Jiang, “Landslide reliability analysis based on transfer coefficient method: A case study from Three Gorges Reservoir,” Journal of Earth Science, vol. 23, no. 2, 187198, 2012. DOI: 10.16285/j.rsm.2006.08.024.
[11] X. W. Hu, H. M. Tang and Y. R. Liu, “Verification of transfer coefficient method applied to landslide stability analysis by physical model test,” Rock and Soil Mechanics, vol. 26, no. 1, 6366, 2005. DOI: 10.16285/j.rsm.2005.01.014.
[12] T. Yang, D. P. Zhou and Y. M. Luo, “A new method considering interlayer interactions for analysis of multilayers landslides,” Chinese Journal of Rock Mechanics and Engineering, vol. 24, no. 7, 11291133, 2005.
[13] X. Huang, C. Wang, T. Wang and Z. Zhang, “Quantification of geological strength index based on discontinuity volume density of rock masses,” International Journal of Heat & Technology, vol. 33, no. 4, 255261, 2015. DOI: 10.18280/ijht.330434.