Analysis of Formation Mechanism of Slightly Inclined Bedding Mudstone Landslide in Coal Mining Subsidence Area Based on Finite–Discrete Element Method

: In this paper, the formation mechanism of a slightly inclined bedding mudstone landslide in the overlying mountain of the coal mining subsidence area of the Tanshan Coal Mine in Ningxia, China, is studied. By means of geotechnical investigation, indoor geotechnical tests, theoretical analysis and other technical means, we ﬁnd the geological environment background of the study area and obtain the physical and mechanical property indexes of the mining landslide in the Tanshan Coal Mine. By combining the numerical simulation of discrete elements and ﬁnite elements, the macro deformation and failure law of the mining mudstone landslide and the displacement and stress nephogram of the failure process are discussed. The results show that the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine is 850 m long from east to west, 500 m wide from north to south and 10,875,000 m 3 in volume. It is composed of Jurassic mudstone and is a traction landslide caused by the coal mining subsidence area. The formation of the landslide is affected by internal factors and inducing factors. The internal factors are mainly geotechnical types and engineering geological properties, and the inducing factors are mainly coal mining activities and rainfall. By analyzing and summarizing the calculation process of the slope model prior to the landslide in 2D-Block and GeoStudio numerical simulation software, the sliding process of the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine is divided into four stages: slope creep, slope deformation, landslide movement and landslide accumulation. GeoStudio software is used to calculate the stability of the Tanshan Coal Mine landslide under natural and rainfall conditions. The landslide is in a stable state under natural conditions and is basically stable under rainfall conditions. By comparing the calculation results of the limit equilibrium method and the ﬁnite element limit equilibrium method, we ﬁnd that the calculated stability coefﬁcient is more accurate when the appropriate constitutive model is selected. The research results have important reference signiﬁcance for the prevention and control of the gently inclined bedding mudstone landslide of the overlying mountain in the coal mining subsidence area of the Loess Plateau.


Introduction
China is rich in coal resources.In recent years, the large-scale development and utilization of coal has brought huge economic and social benefits, providing a basic guarantee for economic development [1,2].Engineering geological disasters caused by coal mining subsidence mainly include geological structure activation; surface cracks, collapse pits and karst collapse; landslides and collapse; rock bursts; mine earthquakes; coal and gas outbursts and other disasters.Mining subsidence is common in the mining area, but when landslide conditions are adequate in the mining subsidence area, mining subsidence in mountainous areas may induce a mining landslide; that is, it will cause the "resurrection" of ancient landslides and the occurrence of new landslides [3].At the same time, the underground coal mining subsidence area weakens the underlying support of the mountain slope, thus inducing the mountain with landslide topographic and geological conditions to produce mining landslide [4,5].For example, the Langshitan landslide caused by the Zijiang coal mine had a length of 2000 m, a width of 5-100 m and large-scale ground fracture zone damage with a visible depth of 5-12 m [3].Generally, after the collapse of a coal mining subsidence area, there are often large cracks, pits and failure forms such as sliding or landslides accompanying the subsidence surface.These characteristics are the obvious signs of landslides in coal mining subsidence areas [6].Coal mining subsidence area landslides not only directly damage the underground and surface construction facilities, resulting in casualties and huge property losses, but also lead to other secondary geological disasters and greater losses [7].Tongchuan Mining Area, with the longest mining time in Shaanxi, has many ground fissures: there are 154 new and old landslides and 351 collapses and falling bodies [3] in the urban area that were caused by mining; the Madaling landslide in Fuxi village, Jiangzhou Town, Duyun City, Guizhou Province, has an overall sliding volume of 190 × 10 4 m 3 due to coal seam mining [8].In addition, unreasonable mining also worsens environmental problems [3].This has brought direct harm to the surrounding ecological environment and the safety of people's lives and property, and it has seriously affected sustainable economic development and social stability in the region [9].Therefore, the study of landslides in coal mining subsidence areas has important guiding significance for disaster prevention and reduction [10].
So far, scholars at home and abroad have conducted extensive and in-depth research on mining landslides and have achieved considerable research results.Mining landslides have mainly been studied from the following three aspects: For the inducing factors of mininginduced landslide, Yang et al. [11] combined the spatiotemporal statistics of landslides and fractures and studied the inducing mechanism of surface subsidence on shallow landslides; Benko and Stead [12], Yin et al. [13] and Xu et al. [14] combined numerical simulation to study the mechanism of landslides caused by changes to the stress field and displacement field in the slope under the influence of goafs.For the deformation characteristics of mining landslides, Long et al. [15] put forward that the failure mode corresponding to the landslide failure stage is the end of goaf subsidence-stress relaxation of rock mass in the slope-creep of soft rock developing creep to the top of the slope-extrusion deformation and failure of the slope; Dong et al. [16] summarized the formation mechanism of landslides as "tensile fracture-creep-shear slip" with a mechanical model through monitoring data analysis and numerical simulation.In addition, for the stability of mining landslides, Long et al. [17] analyzed and summarized the research status of mining landslides and summarized the existing stability analysis methods of mining landslides in combination with theoretical research and numerical analysis; Zhang et al. [18], Salmi et al. [19] and Dai et al. [20] analyzed the influence of mining on slope stability through numerical simulation and limit equilibrium.Yang et al. [11] used interferometric synthetic aperture radar (InSAR) technology and temporal and spatial statistics of landslides and fractures to study the expansion process of slow settlement caused by underground mining and its impact on the occurrence of shallow landslides.Yang et al. [21] put forward a fuzzy comprehensive evaluation model for slope stability in coal mining subsidence areas and, based on this model, evaluated the slope stability in the subsidence zone of the Antaibao Open Pit Coal Mine.Based on the above research, it is found that although adverse geological structural conditions are the internal factors for landslide formation, large-scale underground coal mining is the main reason for the occurrence of mining-induced landslides [7].
So far, based on the theoretical research of traditional landslides, mining subsidence and rock mass mechanics, many methods such as physical model tests, field monitoring and numerical simulation have been used to study landslides in coal mining subsidence areas [15,22,23].As an inversion method, physical model tests are more and more favored by geologists because they can better simulate the actual situation in the field, so as to obtain more accurate evolution information [24].For example, Zhao et al. [25] and Fan et al. [26] analyzed the deformation and failure mechanism of mining landslides through physical model tests; Tang et al. [23] and Dai et al. [27] discussed the deformation and failure characteristics of mining landslides in combination with physical model tests.However, for physical model tests, it is very difficult to design a model test that fully reflects the real prototype because it has high requirements for model framework and materials, needs to consume more human and financial resources, and has a longer research cycle with poor flexibility and convenience [28,29].Tao et al. [30] and Liang and Tang [31] conducted real-time monitoring on the evolution characteristics of landslides based on on-site monitoring, but compared with other research methods, on-site monitoring requires a lot of time, manpower and expensive monitoring instruments [29,32].At present, most researchers study the formation mechanism of landslides based on numerical simulation; for example, Zhao et al. [33], Wang et al. [8] and Salmi et al. [19] have studied the deformation development process of landslides in goaf in combination with numerical simulations and have analyzed their formation mechanisms.Mao et al. [34] analyzed the shallow damage effect of loess slope under the action of rainfall and plant roots by numerical simulation method.Compared with the other two methods, numerical simulations have certain applicability and reliability for the study of coal mining subsidence area landslides, it is not necessary to spend a lot of time and work in the research process, and they are more economical and applicable.Discrete element numerical simulations mainly emphasize the discontinuity of rock masses, allowing the blocks to slip or detach after being pulled, or even free fall from the parent body [35].Finite element numerical simulation follows the continuity assumption of the medium.It replaces an actual structure or continuum with an approximate equivalent physical model composed of multiple interconnected elements.Based on the basic principles of continuum mechanics and the physical characteristics of the elements, a set of equations representing the relationship between force and displacement is established, and the stress, strain and other auxiliary quantities of each element are obtained accordingly [35,36], Therefore, it is not suitable for simulating large deformation and damage fracture phenomena [36].The distinct advantage of discrete element methods is that they can easily handle the mechanical problems of discontinuous media.However, because discrete element methods are not good at describing the plastic deformation of materials, and considering the calculation scale and efficiency, they need to be combined with finite element methods [36].This paper combines the finite element method with the discrete element method to realize the numerical simulation of the formation mechanism of a landslide in a coal mining subsidence area.Li et al. [37] adopted a continuous-discontinuous element method to simulate the in situ fluidized mining filling process of coal.Antolini et al. [38] combined the characteristics of a finite element method and a discrete element method to describe the FDEM simulation results of the trigger mechanism and evolution scenario of the Torgiovannetto Di Assisi rock landslide (in central Italy).
The origin of mudstone is mainly continental sedimentation belonging to soft rock with weak water permeability, weak weathering resistance and easy disintegration in water [39].Different from the bedding failure mode of hard rock landslides, the dip angle of mudstone landslides is usually relatively slow (generally less than 20 • ).According to the Mohr Coulomb failure criterion, it is difficult to start such gently inclined stratum landslides, and it is difficult to form huge sliding thrust and displacement in a short time only depending on the gravity component of the rock mass itself.Many experts and scholars have analyzed the failure mechanism, numerical simulation and other aspects of gently inclined bedding mudstone landslides in combination with engineering practice, and they agree that excavation unloading, weak interlayer zone [40,41], rainfall seepage [42,43] and other factors are the main factors leading to the failure of mudstone slopes [44].The Tanshan Township, Yuanzhou District, Ningxia, China, is rich in mineral resources.Underground mining has induced a large number of mining mudstone landslides in this area, which seriously endangers local economic construction and the safety of people's lives and property.Therefore, taking a slightly inclined bedding mudstone landslides in the Tanshan Coal Mine landslide as the research object, based on geotechnical investigation data and indoor geotechnical test data combined with discrete element numerical simulation and finite element numerical simulation methods, this paper systematically studies the formation mechanism and stability of the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine.It is not only of great significance to the prevention and control of geological disasters in Yuanzhou District, Guyuan City, Ningxia, China, but also provides reference materials for the study of similar landslides.

Regional Geological Environment Background
The Tanshan Coal Mine landslide in Ningxia, China is located in Group 2, Tanshan Village, Tanshan Township, Yuanzhou District, Guyuan City, about 1 km away from the Tanshan Township government.The coordinates are E 106 • 18 18.60 and N 36 • 23 8.40 .The area where the Tanshan Coal Mine landslide is located is a typical continental semi-arid climate, characterized by less annual precipitation, concentrated rainy season, uneven precipitation distribution, strong evaporation, low annual average temperature, large temperature difference between day and night and more disastrous weather.In addition, the area is poor in water resources and poor in water quality.The groundwater type is mainly stagnant water in the upper layer, which is mostly between mudstone and loess, and the water level and water volume change obviously with seasons.The landslide area of the Tanshan Coal Mine is a loess hilly area with large topographic relief.Most of them are based on Jurassic or tertiary red mudstone.The Quaternary loess and loess-like soil fluctuate along and cover the ancient terrain with a thickness of less than 50 m.The terrain cutting is relatively strong, with an altitude of 1600-2000 m.It is dominated by loess beams and coexisting loess beams and loess hills, and the gully is in a "V" shape.The lower mudstone is often exposed at the bottom of the valley.The geotectonic structure in the area is located in the North Qilian fold area.At the end of Neogene, the tectonic line in the area extends in a north-northwest direction due to strong influence by the West Malaya movement and strong compression from the southwest.The mine was built in 1956.It is a small state-owned enterprise with a designed production capacity of 90,000 tons yearly and an actual production capacity of 107,000 tons yearly.The mining method is underground mining.The landslide is located in the west of Huangtuliang, Group 2, Tanshan Village.The front edge is a deep gully, which is nearly north-south, is high in the south and low in the north, and has strong cutting.The back wall of the landslide is close to the top of the beam.The terrain is relatively narrow, and the sliding mass is gentle.It is generally divided into two platforms, which are steep at the top and gentle at the bottom.The overall slope is about 10 • , and the sliding direction is 290 • .The original slope height was about 180 m, and the slope was about 11 • .The landslide morphology of the Tanshan Coal Mine is shown in Figure 1.The geological map and profile of the landslide in the coal mining subsidence area of Tanshan, China, are shown in Figure 2.

Laboratory Geotechnical Test
In order to obtain numerical simulation calculation parameters, density testing and direct shear testing were carried out according to the Chinese "Standard for Soil Test Method" (GB/T 50123-2019).The landslide soil samples collected by engineering geological drilling were carefully placed and brought back to the laboratory.The undisturbed soil samples were prepared with a ring knife with a diameter of 61.8 mm and a height of 20.0 mm.
The density test adopted the ring knife method, and the density of the landslide soil sample was obtained through Formula (1).
where ρ is the sample density (g/cm 3 ), m 0 is the mass of the soil sample (g), and V is the ring-cutter volume (cm 3 ).
The direct shear test adopted the ZJ strain-controlled direct shear instrument produced by the Nanjing Soil Instrument Factory in China for rapid shear tests.The test is designed to be carried out under vertical pressures of 50 kPa, 100 kPa, 200 kPa and 300 kPa, and the tests were repeated twice under each vertical pressure.During the shear process, the shear rate was set as 0.8 mm/min, and the test was stopped when the shear deformation reached 6 mm.We calculated the shear stress of the soil sample and drew the relationship curve between shear stress τ and shear displacement ∆L, shear strength S and the vertical unit pressure p, obtained shear strength S, cohesion c and the internal friction angle ϕ [45].

Discrete Element Numerical Simulation Method (1) Model establishment
In order to analyze the key problems studied in this paper and avoid the lengthy calculation process, we make appropriate simplified assumptions before the simulation: 1 Because the scale of the slope body before the landslide in the Tanshan Coal Mine is much larger than the scale of the fluctuation of the slope shape, the original slope shape is simplified to a straight slope shape. 2 Because the interface between moderately weathered mudstone and slightly weathered mudstone is difficult to determine, and the difference of engineering properties is small, it has little impact on the problems to be studied.Therefore, all the rock layers under the strongly weathered mudstone layer are collectively called a moderately weathered mudstone layer.
According to the actual situation of a slightly inclined bedding mudstone landslide in the Tanshan Coal Mine landslide, the main section of the landslide is selected to establish the slope discrete element model.The model is 2000 m wide and 800 m high, as shown in Figure 3.The boundary conditions of the slope discrete element method are that the left, right and lower boundaries are fixed boundaries in the X-and Y-directions, and other boundaries are free boundaries, as shown in Figure 3.
(2) Constitutive model The 2D-Block software provides two constitutive models: angle-edge contact and a combination of angle-edge and edge-edge contact.In this paper, the constitutive models of combined angle-edge and edge-edge contact are selected.The model is an elastic-plastic tension-free contact model that is speed-independent.We assumed no tension between blocks, and that plastic shear slip will occur when the tangential force F s reaches a certain maximum, which is determined by the following formula: where F n is the normal force, ϕ is the internal friction angle, (F s ) max is the maximum tangential force and C is cohesion.(2) Constitutive model The 2D-Block software provides two constitutive models: angle-edge contact and a combination of angle-edge and edge-edge contact.In this paper, the constitutive models of combined angle-edge and edge-edge contact are selected.The model is an elastic-plastic tension-free contact model that is speed-independent.We assumed no tension between blocks, and that plastic shear slip will occur when the tangential force s F reaches a cer- tain maximum, which is determined by the following formula: where n F is the normal force, ϕ is the internal friction angle, ( ) max s F is the maximum tangential force and c is cohesion.Assuming that the length and width of two contact blocks are a and b, respec- tively, the elastic modulus is E and the Poisson's ratio is γ, according to the theory of elasticity: Thus, the normal stiffness coefficient is obtained: The tangential stiffness coefficient s K can be obtained from the normal stiffness co- For angle-edge contact, the above formula can be directly used for calculation.In this constitutive model, edge-edge contact exists as an independent model.When the contact side length of two blocks is L, the uniformly distributed contact stiffness on the side is n JK and s JK , and the uniformly distributed bonding force is Jc , the total stiffness on the side is: Assuming that the length and width of two contact blocks are a and b, respectively, the elastic modulus is E and the Poisson's ratio is γ, according to the theory of elasticity: Thus, the normal stiffness coefficient is obtained: The tangential stiffness coefficient K s can be obtained from the normal stiffness coefficient K n : For angle-edge contact, the above formula can be directly used for calculation.In this constitutive model, edge-edge contact exists as an independent model.When the contact side length of two blocks is L, the uniformly distributed contact stiffness on the side is JK n and JK s , and the uniformly distributed bonding force is Jc, the total stiffness on the side is: For edge-edge contact, a set of joint constants JK n , JK s , Jc and J ϕ shall be given during calculation, which shall be transformed into the total stiffness on the edge according to Formulas (6)-( 8) and comply with the sliding criterion (Formula (2)).When the joint parameter JK n of edge-edge contact and the joint parameter K n of edge-angle contact meet JK n • L < K n , the edge-edge contact will be automatically transformed into angle-edge contact.
(3) Parameter selection When using 2D-Block discrete element software to simulate slope deformation, it is necessary to provide the stiffness value and strength index of the rock and soil mass.In this paper, the model is established according to the actual situation of the geological section, and the discrete elements are divided according to the accuracy requirements of numerical simulation.The Loess and strongly weathered mudstone are divided into elements according to side length a = 10 and b = 10; the moderately weathered mudstone is divided into elements according to side length a = 20 and b = 20; and the elastic modulus of loess is E = 6 × 10 5 Pa.Substituting into Formulas (4)- (7), we obtain K n = (Ea)/(2b) = 3 × 10 5 Pa, 4 Pa, and other soil layer parameters can be obtained in the same way.Due to the discretization of the slope, the cohesion c is taken as 0 and ϕ is taken as the internal friction angle corresponding to the residual strength.The selection of stiffness parameters and strength indexes is shown in Table 1.Firstly, we establish a goaf-free model.Then, according to the geological profile of a slightly inclined bedding mudstone landslide in the Tanshan Coal Mine after mining, a calculation model (as shown in Figure 4) is established for the critical state before sliding, that is, the state when the through weak zone is formed but there was no severe sliding, and it is divided into elements.The parts with complex stress changes that require key analysis are densified by 0.5 times and divided into a total of 11,420 nodes and 11,335 elements.For different models, there is no fixed grid density for reference, but the rationality of grid density can be determined by a grid convergence test.The boundary conditions in Figure 4 are that the left and right boundaries are X-direction displacement constraints, the lower boundary is X-and Y-direction displacement constraints, and the upper boundary is a free boundary.The model is taken as an approximate ideal elastic-plastic model, and the Mohr-Coulomb yield criterion is selected as the failure criterion.(2) Parameter selection The calculation parameters of landslide finite element numerical simulation are the parameters of each rock stratum of the slope when the landslide does not occur, as shown in Table 2: (

2) Parameter selection
The calculation parameters of landslide finite element numerical simulation are the parameters of each rock stratum of the slope when the landslide does not occur, as shown in Table 2: The purpose of landslide stability analysis is to judge the stability state of a landslide and provide stability analysis data for landslide treatment.The stability degree of a landslide is usually expressed by stability coefficient F s [46,47].At present, the most widely used methods are limit equilibrium and numerical analysis.
The limit equilibrium method assumes that the landslide is a rigid body under the condition of a certain instability model.Through study of the stress state and effective strength of the instability separation surface, the stability coefficient of the landslide is calculated according to the limit equilibrium principle to characterize the stability of the landslide [48].At present, many different methods have been developed from the limit equilibrium method, and all of them are basically similar.The difference is that different equilibrium equations are included and satisfied and they have different assumptions on the relationship between inter-slice forces and conditional forces, as shown in Table 3.In the actual engineering analysis, when the limit equilibrium method is adopted, the method that meets the balance of force and moment should be adopted as much as possible, and some safety factor equations with stricter mathematical significance can be easily solved using computer software.As for numerical analysis methods, they mainly use some method to calculate the stress distribution and deformation of a slope or landslide, study the change process of stress and strain in rock and soil masses, and obtain the local stability coefficient at each point to judge the stability of a slope or landslide.For stability analysis by a finite element method, firstly, the stress distribution in the site needs to be obtained by finite element theoretical analysis, and then, the calculated finite element stress is introduced into the traditional limit equilibrium equation (China Imitation Technology Co 2011).The stability coefficient can be calculated by calculating the sliding force and anti-sliding force of each soil strip.The formula of stability coefficient F s is defined as: where S r is the total shear resistance along the whole sliding surface and S m is the total sliding force along the whole sliding surface.

Limit Equilibrium Method and Finite Element Limit Equilibrium Method Based on GeoStudio
In this paper, the limit equilibrium method of the SLOPE/W module of GeoStudio software and the finite element limit equilibrium method of the coupling function of SIGMA/W and SLOPE/W modules are used to calculate the post-sliding stability of the Tanshan Coal Mine under natural and rainfall conditions, and the stability coefficients calculated by different methods are compared and analyzed.The calculation model and parameter selection of the limit equilibrium method and the finite element limit equilibrium method are as follows: (1) Model establishment According to the geological profile obtained from field investigation, the calculation model of the limit equilibrium method is established.In the finite element analysis, 4673 elements and 4819 nodes are divided according to the program.The constraint conditions are that the left and right boundaries are X-direction displacement constraints, the lower boundary is X-and Y-direction displacement constraints, and the upper boundary is a free boundary.
(2) Parameter selection Under natural working conditions, the natural unit weight and shear strength indexes obtained from direct shear testing are taken for each layer.Under rainfall conditions, due to the small permeability coefficient of loess, it is proposed to take a rainfall infiltration depth of 2 m, which takes the saturated unit weight and the shear strength index obtained from the saturated direct shear test; at the same time, the shear strength indexes obtained from the saturated unit weight and the saturated direct shear test are also taken for the mudstone in the slip zone (as shown in Table 4).According to the geotechnical engineering investigation, the length of the landslide in the east-west direction is about 850 m, the maximum width to the north and south is 500 m, the thickness of the landslide is 10-20 m, the plane area is 475,000 m 2 , the volume is 10,875,000 m 3 , and the sliding direction is 295 • .From the material composition and scale of the landslide, it is a giant bedrock landslide.From analysis of the landslide deformation characteristics and motion properties, it is a traction landslide.The main body of the landslide in the Tanshan Coal Mine is composed of Jurassic mudstone.Reconstructed by history and human activities reveals that the original topography of the landslide has been destroyed and changed, but the overall shape of the landslide is still basically clear, the outline of the landslide can be identified, the shape of the back wall of the landslide is irregular, the height of the scarp is 1-5 m, and the boundary line of the side wall is obvious.There are two levels of platforms on the sliding body, which are steep up and slow down, with an overall slope of about 20 • .The leading edge of the landslide was in the air and was a deep-cut valley.The original slope height was about 180 m, and the slope was about 30 • .

Deformation and Failure Characteristics of the Landslide
The mining point of the Tanshan Coal Mine is located in front of the landslide.Due to the influence of underground mining, the front edge of the landslide subsidence occurs, making the whole slope a two-stage platform.The front platform has been used to build factories and for other operations of the Tanshan Coal Mine after being flattened, and the rear platform has been mostly used to build houses or to develop farmland.
The colliery landslide of the Tanshan Coal Mine occurred in 1966; in recent years, because of rainfall, coal mining subsidence, earthquakes and other comprehensive factors, there have been signs of revival.At present, there is a set of tensile cracks in the trailing edge of the landslide that strike 32 • , and the dominant crack is about 300 m long.Most of the cracks are buried by the local residents, and some are still identifiable.There are several small water holes that are 0.5-2 m deep.Rainwater flows from cracks during rainfall, which makes the back of the landslide unstable.Due to the continuous erosion of the gully in front of the landslide, a steep free face is formed.
Affected by sliding, the original structure of the sliding body's rock and earth mass is disturbed, the structure is loose, and cracks and holes develope.The sliding body is scoured by surface water, forming a gully of different sizes, and the leading edge has been scoured by water in the gully for a long time, forming a free face.There are tensile cracks in the trailing edge.In addition, affected by continuous rainfall in the rainy season, the moisture content of the slope increases continuously.On the one hand, this increases the weight of the slope, and on the other hand, rainwater infiltration reduces the shear strength of the sliding zone.Under the combined action of these factors, the landslide has had a trend of revival in recent years; this must be given adequate attention.

Numerical Simulation Results of Discrete Element Method
(1) Numerical simulation results only considering the coal mining subsidence area The discrete element numerical simulation under the action of the coal mining subsidence area has 283,401 iteration steps.Figure 5 shows that the stress state of the rock mass changes under the action of the coal mining subsidence area, and the rock strata above the two ends of the coal mining subsidence area are compressed, which makes the support pressure transfer to the depth of both sides.The roof and floor of the working face in the coal mining subsidence area form a decompression area, and the lower strata have small elastic recovery deformation.Under the action of gravity, the roof strata bends, breaks and falls in the coal mining subsidence area direction and forms a layered subsidence, and there are obvious delamination, crack and fracture phenomena.Displacement is transmitted to the surface, forming a subsidence basin on the surface, and coal mining subsidence area cracks appear at the edge of the subsidence basin and at the places with low rock strength and soil mass.These factors can accelerate the creep failure of the slope.
Figure 5 shows that there are obvious tensile cracks in the trailing edge of the slope, and the maximum displacement of the block is 1254 m.That is to say, even if the support fails completely and the coal mining subsidence area collapses completely, the maximum displacement of the trailing edge of the slope is only 1254 m, which is far from the critical state that can cause the slope slide and is inconsistent with the field survey.It can be judged that the influence of the coal mining subsidence area alone is not enough to form such a large landslide on a gentle slope.
breaks and falls in the coal mining subsidence area direction and forms a layered subsidence, and there are obvious delamination, crack and fracture phenomena.Displacement is transmitted to the surface, forming a subsidence basin on the surface, and coal mining subsidence area cracks appear at the edge of the subsidence basin and at the places with low rock strength and soil mass.These factors can accelerate the creep failure of the slope.Figure 5 shows that there are obvious tensile cracks in the trailing edge of the slope, and the maximum displacement of the block is 1254 m.That is to say, even if the support fails completely and the coal mining subsidence area collapses completely, the maximum displacement of the trailing edge of the slope is only 1254 m, which is far from the critical state that can cause the slope slide and is inconsistent with the field survey.It can be judged that the influence of the coal mining subsidence area alone is not enough to form such a large landslide on a gentle slope.
(2) Numerical simulation results of landslide process The slightly inclined bedding mudstone landslide in the Tanshan Coal Mine was formed in the rainy season of 1966, and underground mining destroyed the stress field of the slope body, resulting in tensile cracks at the rear edge of the slope body and collapse basins and mining cracks at the front edge.Rainfall and surface water infiltrated downward along the rock mass cracks formed by tensile cracks, mining cracks and mining subsidence.Since the strength and permeability of strongly weathered mudstone are smaller than those of moderately weathered mudstone, water easily combined with mudstone at the interface between moderately weathered mudstone and strongly weathered mudstone after rainfall and surface water infiltration, resulting in saturated softening of mudstone, reducing the shear strength of rock and soil and forming local weak zones.At the (2) Numerical simulation results of landslide process The slightly inclined bedding mudstone landslide in the Tanshan Coal Mine was formed in the rainy season of 1966, and underground mining destroyed the stress field of the slope body, resulting in tensile cracks at the rear edge of the slope body and collapse basins and mining cracks at the front edge.Rainfall and surface water infiltrated downward along the rock mass cracks formed by tensile cracks, mining cracks and mining subsidence.Since the strength and permeability of strongly weathered mudstone are smaller than those of moderately weathered mudstone, water easily combined with mudstone at the interface between moderately weathered mudstone and strongly weathered mudstone after rainfall and surface water infiltration, resulting in saturated softening of mudstone, reducing the shear strength of rock and soil and forming local weak zones.At the same time, due to the continuous collapse of the coal mining subsidence area, the cracks in the rock and soil mass further increased, and the process continued to repeat.Local softening eventually formed a penetrating weak zone.At the same time, the weight of the upper rock and soil mass increase after water absorption.The combined action of these factors eventually led to the sliding failure of the slope.
On the basis of the original model, a 1 m thick weak zone is set on the interface between the strongly weathered and moderately weathered mudstone to simulate the sliding process of the landslide.Self-adaptive damping is used in the solution control, and the iterative step length is 0.000146 s.Through iterative calculation, the simulation starts from motion to stability (iterative completed), and the total number of iterations is 607,693.The whole process of slightly inclined bedding mudstone landslide sliding in the Tanshan Coal Mine is basically restored.The simulation process is shown in Figures 6-9.
When the numerical simulation iteration reaches 10,000 steps (Figure 6), the rock mass above the coal mining subsidence area deforms, and the deformation is transmitted to the surface.A small collapse can be seen on the surface, and mining cracks are formed on the contact surface of strong weathering and medium weathering, and the rock and soil mass is loose.Due to the influence of the frontal collapse basin, the slope begins to decline, and tensile cracks appear at the rear edge.The maximum displacement of the rear edge block was 0.5831 m according to the block tracking function of 2D-Block, and the displacement gradually increased from bottom to top in the vertical direction of the slope.
tween the strongly weathered and moderately weathered mudstone to simulate the sliding process of the landslide.Self-adaptive damping is used in the solution control, and the iterative step length is 0.000146 s.Through iterative calculation, the simulation starts from motion to stability (iterative completed), and the total number of iterations is 607,693.The whole process of slightly inclined bedding mudstone landslide sliding in the Tanshan Coal Mine is basically restored.The simulation process is shown in Figures 6-9.When the numerical simulation iteration reaches 10,000 steps (Figure 6), the rock mass above the coal mining subsidence area deforms, and the deformation is transmitted to the surface.A small collapse can be seen on the surface, and mining cracks are formed on the contact surface of strong weathering and medium weathering, and the rock and soil mass is loose.Due to the influence of the frontal collapse basin, the slope begins to decline, and tensile cracks appear at the rear edge.The maximum displacement of the rear edge block was 0.5831 m according to the block tracking function of 2D-Block, and the displacement gradually increased from bottom to top in the vertical direction of the slope.
When the numerical simulation iteration is 10,000 to 200,000 steps (Figure 7), the coal mining subsidence area completely collapses, and the displacement of the rock layer is gradually transferred to the surface, forming a subsidence basin on the surface.The mining cracks caused by the coal mining subsidence area and the movement of the slope formed tensile cracks at the trailing edge.The upper rock mass of the slope slides under the action of self-weight stress, and the tensile cracks at the trailing edge are further widened and developed to reach the deep areas.The sliding body and the trailing edge are gradually separated, and the vertical displacement and horizontal displacement increase.Small-scale collapses have occurred in the front face.x FOR PEER REVIEW 14 of 27 During the iteration of 300,000-450,000 steps, after a long-term creep deformation stage, the coal mining subsidence area collapse makes the supporting effect of the lower rock mass on the upper rock mass gradually decrease, the tensile crack at the rear edge increases continuously, and the maximum displacement reaches 15.231 m.The sliding body slides violently along the sliding surface, shears out along the shear outlet, and slides towards the leading edge gully.The maximum displacement of the block falling into the gully at the leading edge reaches 20.334 m and finally forms a landslide (Figure 8).During the iteration of 300,000-450,000 steps, after a long-term creep deformation stage, the coal mining subsidence area collapse makes the supporting effect of the lower rock mass on the upper rock mass gradually decrease, the tensile crack at the rear edge increases continuously, and the maximum displacement reaches 15.231 m.The sliding body slides violently along the sliding surface, shears out along the shear outlet, and slides towards the leading edge gully.The maximum displacement of the block falling into the gully at the leading edge reaches 20.334 m and finally forms a landslide (Figure 8).For the iteration of 450,000-607,693 steps, the landslide, after severe sliding along with the release of energy, slides slowly and gradually, and the landslide body accumulates in front of the landslide and the slope gully (Figure 9).

Finite Element Numerical Simulation Results
(1) Displacement analysis In Figure 10, the green part is the model grid, and the size and direction of the red arrow represents the displacement vector of the grid node.It can be seen from Figure 10 that the displacement of the slope mainly occurs in the middle and lower part of the slope and above the coal mining subsidence area.Due to the influence of mining, the stress field of the rock mass above the coal mining subsidence area is changed, and the stress of the rock mass in the mining area is released.The rock mass moves downward under selfweight stress and the extrusion of overlying strata, and the coal mining subsidence area is filled.Displacement is transmitted to the surface, and subsidence basins are formed on the surface.It can also be seen from the figure that the weak zone plays an important role in the formation of the landslide.The displacement vector has an obviously offset dislocation above and below the weak zone, and the displacement of rock soil above the weak zone is significantly greater than that below.This shows that under the influence of mining, the upper and lower parts of the rock and soil mass in the weak zone have downward displacement.Due to the effect of the weak zone, the displacement of the upper rock and soil mass is significantly greater than that of the lower rock and soil mass, and this difference in displacement makes the upper rock and soil mass slide in a relative manner.When the numerical simulation iteration is 10,000 to 200,000 steps (Figure 7), the coal mining subsidence area completely collapses, and the displacement of the rock layer is gradually transferred to the surface, forming a subsidence basin on the surface.The mining cracks caused by the coal mining subsidence area and the movement of the slope formed tensile cracks at the trailing edge.The upper rock mass of the slope slides under the action of self-weight stress, and the tensile cracks at the trailing edge are further widened and developed to reach the deep areas.The sliding body and the trailing edge are gradually separated, and the vertical displacement and horizontal displacement increase.Small-scale collapses have occurred in the front face.
During the iteration of 300,000-450,000 steps, after a long-term creep deformation stage, the coal mining subsidence area collapse makes the supporting effect of the lower rock mass on the upper rock mass gradually decrease, the tensile crack at the rear edge increases continuously, and the maximum displacement reaches 15.231 m.The sliding body slides violently along the sliding surface, shears out along the shear outlet, and slides towards the leading edge gully.The maximum displacement of the block falling into the gully at the leading edge reaches 20.334 m and finally forms a landslide (Figure 8).
For the iteration of 450,000-607,693 steps, the landslide, after severe sliding along with the release of energy, slides slowly and gradually, and the landslide body accumulates in front of the landslide and the slope gully (Figure 9).

Finite Element Numerical Simulation Results
(1) Displacement analysis In Figure 10, the green part is the model grid, and the size and direction of the red arrow represents the displacement vector of the grid node.It can be seen from Figure 10 that the displacement of the slope mainly occurs in the middle and lower part of the slope and above the coal mining subsidence area.Due to the influence of mining, the stress field of the rock mass above the coal mining subsidence area is changed, and the stress of the rock mass in the mining area is released.The rock mass moves downward under self-weight stress and the extrusion of overlying strata, and the coal mining subsidence area is filled.Displacement is transmitted to the surface, and subsidence basins are formed on the surface.It can also be seen from the figure that the weak zone plays an important role in the formation of the landslide.The displacement vector has an obviously offset dislocation above and below the weak zone, and the displacement of rock soil above the weak zone is significantly greater than that below.This shows that under the influence of mining, the upper and lower parts of the rock and soil mass in the weak zone have downward displacement.Due to the effect of the weak zone, the displacement of the upper rock and soil mass is significantly greater than that of the lower rock and soil mass, and this difference in displacement makes the upper rock and soil mass slide in a relative manner.
weight stress and the extrusion of overlying strata, and the coal mining subsidence area is filled.Displacement is transmitted to the surface, and subsidence basins are formed on the surface.It can also be seen from the figure that the weak zone plays an important role in the formation of the landslide.The displacement vector has an obviously offset dislocation above and below the weak zone, and the displacement of rock soil above the weak zone is significantly greater than that below.This shows that under the influence of mining, the upper and lower parts of the rock and soil mass in the weak zone have downward displacement.Due to the effect of the weak zone, the displacement of the upper rock and soil mass is significantly greater than that of the lower rock and soil mass, and this difference in displacement makes the upper rock and soil mass slide in a relative manner.It can be seen from Figure 11 that the displacement of the slope in the X-direction under the influence of the coal mining subsidence area is very different, with the maximum value of 22.3 m and the minimum value of −5.2 m.The lower part of the slope is It can be seen from Figure 11 that the displacement of the slope in the X-direction under the influence of the coal mining subsidence area is very different, with the maximum value of 22.3 m and the minimum value of −5.2 m.The lower part of the slope is mainly due to the influence of the coal mining subsidence area, and the direction of the displacement in the lower part of the subsidence basin points to the coal mining subsidence area, so the X-direction displacement is negative.For the upper slope under the influence of the coal mining subsidence area, the displacement direction points to the coal mining subsidence area, so the X-direction displacement value is positive.At the same time, due to the relationship with the weak zone, the X-direction of the rock mass in the upper and lower parts of the weak zone has obvious dislocation, which shows that the displacement values in the X-direction of the rock mass in the upper and lower parts of the weak zone are very different.The displacement of the upper rock mass along the weak zone to the slope toe is much larger than that of the lower rock mass.mainly due to the influence of the coal mining subsidence area, and the direction of the displacement in the lower part of the subsidence basin points to the coal mining subsidence area, so the X-direction displacement is negative.For the upper slope under the influence of the coal mining subsidence area, the displacement direction points to the coal mining subsidence area, so the X-direction displacement value is positive.At the same time, due to the relationship with the weak zone, the X-direction of the rock mass in the upper and lower parts of the weak zone has obvious dislocation, which shows that the displacement values in the X-direction of the rock mass in the upper and lower parts of the weak zone are very different.The displacement of the upper rock mass along the weak zone to the slope toe is much larger than that of the lower rock mass.It can be seen from Figure 12 that the Y-direction displacement of the slope is mainly affected by the coal mining subsidence area.Due to stress release, the rock layer above the coal mining subsidence area has downward displacement with a maximum value of 4.8 m, and the displacement value decreases with the increase in elevation.This indicates that due to the strength effect of the rock stratum itself, the displacement of the rock mass caused by the coal mining subsidence area decreases continuously in the process of transferring to the surface.The influence of the weak zone on the Y-direction displacement is much smaller than that of the coal mining subsidence area, but the dislocation of the contour can also be seen in the diagram, which shows that the rock mass in the upper part of the weak zone slides along the weak zone, and that the X-direction displacement is primary while the Y-direction displacement is secondary.It can be seen from Figure 12 that the Y-direction displacement of the slope is mainly affected by the coal mining subsidence area.Due to stress release, the rock layer above the coal mining subsidence area has downward displacement with a maximum value of 4.8 m, and the displacement value decreases with the increase in elevation.This indicates that due to the strength effect of the rock stratum itself, the displacement of the rock mass caused by the coal mining subsidence area decreases continuously in the process of transferring to the surface.The influence of the weak zone on the Y-direction displacement is much smaller than that of the coal mining subsidence area, but the dislocation of the contour can also be seen in the diagram, which shows that the rock mass in the upper part of the weak zone slides along the weak zone, and that the X-direction displacement is primary while the Y-direction displacement is secondary.
It can be seen from Figure 12 that the Y-direction displacement of the slope is mainly affected by the coal mining subsidence area.Due to stress release, the rock layer above the coal mining subsidence area has downward displacement with a maximum value of 4.8 m, and the displacement value decreases with the increase in elevation.This indicates that due to the strength effect of the rock stratum itself, the displacement of the rock mass caused by the coal mining subsidence area decreases continuously in the process of transferring to the surface.The influence of the weak zone on the Y-direction displacement is much smaller than that of the coal mining subsidence area, but the dislocation of the contour can also be seen in the diagram, which shows that the rock mass in the upper part of the weak zone slides along the weak zone, and that the X-direction displacement is primary while the Y-direction displacement is secondary.(

2) Stress analysis
It can be seen that the stress field distribution of the slope before coal mining fully conforms to the distribution law of a self-weight stress field.With the increase in depth, the stress contour distribution is generally consistent with the surface morphology, and the slope stability is very high.
After the failure of the coal mining subsidence area support, the stress field of the slope changed dramatically.Figure 13 shows that in a position far from the coal mining subsidence area, the contour of the major principal stress is still nearly parallel to the slope, while in a position near the coal mining subsidence area or the weak zone, the stress concentration occurs in many places, and the contour is nearly perpendicular to the slope.In the coal mining subsidence area boundary and the upper supporting pressure area, due to compression of the upper rock mass, the stress concentration is obvious, and the maximum value reaches 34.78 MPa.In the fully mining area above the middle of the coal mining subsidence area, the major principal stress changes from compressive stress to tensile stress, forming a tensile stress area, and the minimum value reaches 9.301 MPa.(2) Stress analysis It can be seen that the stress field distribution of the slope before coal mining fully conforms to the distribution law of a self-weight stress field.With the increase in depth, the stress contour distribution is generally consistent with the surface morphology, and the slope stability is very high.
After the failure of the coal mining subsidence area support, the stress field of the slope changed dramatically.Figure 13 shows that in a position far from the coal mining subsidence area, the contour of the major principal stress is still nearly parallel to the slope, while in a position near the coal mining subsidence area or the weak zone, the stress concentration occurs in many places, and the contour is nearly perpendicular to the slope.In the coal mining subsidence area boundary and the upper supporting pressure area, due to compression of the upper rock mass, the stress concentration is obvious, and the maximum value reaches 34.78 MPa.In the fully mining area above the middle of the coal mining subsidence area, the major principal stress changes from compressive stress to tensile stress, forming a tensile stress area, and the minimum value reaches 9.301 MPa.
In the upper part of the slope near the weak zone, the major principal stress appears as wave-like changes, which shows the stress concentration near the weak zone.This stress concentration led to the failure of rock and soil in the weak zone.The negative value of major principal stress appears at the rear edge of the slope, which indicates that under the action of the coal mining subsidence area, the collapse basin and mining-induced cracks appear at the front edge of the slope, leading to downward sliding of the main sliding section.The rear edge of the slope loses its support, which causes tensile failure, resulting in tensile cracks.The range of this area is roughly consistent with the position of tensile cracks in the investigation.

Calculation Results of Landslide Stability
The calculation results of the landslide stability coefficient are shown in Table 5: Table 5. Calculation results of landslide stability coefficient after sliding in the Tanshan Coal Mine.In the upper part of the slope near the weak zone, the major principal stress appears as wave-like changes, which shows the stress concentration near the weak zone.This stress concentration led to the failure of rock and soil in the weak zone.The negative value of major principal stress appears at the rear edge of the slope, which indicates that under the action of the coal mining subsidence area, the collapse basin and mining-induced cracks appear at the front edge of the slope, leading to downward sliding of the main sliding section.The rear edge of the slope loses its support, which causes tensile failure, resulting in tensile cracks.The range of this area is roughly consistent with the position of tensile cracks in the investigation.

Calculation Results of Landslide Stability
The calculation results of the landslide stability coefficient are shown in Table 5: The post-sliding stability coefficients in the natural state calculated by the limit equilibrium method and the finite element method are both greater than 1.15, which belongs to the stable state and has small risk.This is consistent with the actual survey; because the coal mining subsidence area has completely collapsed, the stability of the landslide has lost its influence.At the same time, after sliding, energy is released, the slope becomes smaller, slope deformation is slow and the stress distribution of the slope is more conducive to the stability of the slope.
Under the rainfall condition, the stability coefficient is lower than that under the natural condition, which fluctuates within the range of 1.05-1.15,which belongs to the basically stable state but still has certain risks.This is because the landslide structure is relatively loose after sliding, and there must be cracks in the trailing edge of the landslide.Under rainfall conditions, rainwater can easily infiltrate along the cracks.On the one hand, this increases the self-weight of the rock and soil; on the other hand, it reduces the shear strength of the sliding zone and reduces the stability coefficient of the landslide.

Necessity of Combining Discrete Element Method with Finite Element Method to Analyze Landslide Formation Mechanism
The discrete element method (DEM) mainly emphasizes the discontinuity of rock masses.The problem domain is composed of many rock mass elements, but these elements do not require complete close contact.The elements not only consist of surface contact or surface-to-point contact, but also allow sliding or separation between blocks after tension and even free-fall away from the parent body [35].The finite element method (FEM) is the most mature and widely used numerical calculation method.The theoretical basis of the finite element method is the principle of minimum potential energy, which replaces an actual structure or continuum with an approximate equivalent physical model composed of multiple interconnected elements.Through the basic principles of structure and continuum mechanics and the physical characteristics of elements, equations representing the relationship between force and displacement are established.Basic unknown physical quantities are obtained by solving the equations, and the stress, strain and other auxiliary values of each element are obtained [35,36].The finite element method possesses the advantage that it can be used for heterogeneous problems, nonlinear materials and various types of problems, can calculate stress and deformation, seepage, consolidation, rheology, wetting deformation, dynamic and temperature problems, and can adapt to complex boundary conditions.However, it also has the disadvantages of large data preparation workload, easy errors with original data and inability to ensure the continuity of some physical quantities in the whole region [49], and finite element numerical simulation software is generally not suitable for simulating large deformation and damage fracture phenomena due to the continuity assumption of the medium and reliance on grid-based calculations [36].It is difficult to effectively simulate the discreteness of rock and soil particles, and it is impossible to consider the complex force conduction between particles and particle displacement [50].The outstanding advantage of the discrete element method is that it can easily deal with the mechanical problems of discontinuous media, discretizing discontinuous bodies or continuum into independent "micro elements" or "particles" with certain physical significance connected through virtual springs.Therefore, it has natural advantages in dealing with fracture and fragmentation problems.However, because the discrete element method is not good at describing the plastic deformation of materials, at the same time, it needs to be combined with the finite element method for the sake of calculation scale and efficiency [36].
Many domestic and foreign scholars have studied the combination of finite and discrete element methods to carry out research on a variety of sites or cases, such as earthquake landslides and rainfall landslides.Hung et al. [51] studied the complex coseismic process of the Aso Bridge landslide during the main earthquake of the Kumamoto earthquake in 2016 through finite and discrete element analysis of the vertical seismic acceleration (VSA).Hung et al. [52] used finite element (FEA) and discrete element (DEA) methods to discuss the driving mechanism before landslide failure and the dynamic jump process after landslide failure from the three aspects of pore water pressure, saturation and displacement of the sliding zone.This shows that the coupling of FEA and DEA can be used to study the impact of rainfall infiltration on landslides and provide useful insights.Hamdi et al. [53] used the finite element/discrete element modeling method in combination with the discrete fracture network (FDEM-DFN) to analyze the ground settlement of the Kiruna hanging wall, focusing on the impact of discontinuity persistence and spacing, which improved the authors' understanding of the extent and mechanism of hanging wall damage with the advance of the mine.In this paper, in order to analyze the influencing factors and formation mechanism of the landslide in the Tanshan Coal Mine more efficiently and clearly, the discrete element method and finite element method are combined.The 2D-Block discrete element software is used to calculate the deformation of landslide in the Tanshan Coal Mine and to analyze the macro deformation and failure law of the slope, and the SIGMA/W module in GeoStudio finite element software is used to calculate the displacement and stress nephogram of the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine, and both are used to simulate and analyze the stress and strain of the slope before sliding; the formation mechanism of the slope before sliding is verified by mutual comparison.

Analysis of the Influencing Factors on the Landslide in the Coal Mining Subsidence Area
Coal mines in Tanshan Township, Ningxia, China, are rich in resources and have a long mining history.The landslides in the area where the Tanshan Coal Mine is located have the characteristics of low angle, large scale and deep sliding, and are formed under the joint action of underground mining of coal seams and a variety of comprehensive factors.The internal influencing factors of the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine in Ningxia, China, are geotechnical types and engineering geological properties, and the inducing influencing factors are coal mining activities and rainfall.

Internal Factors
The internal influencing factors of the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine in Ningxia, China, are geotechnical types and engineering geological properties.The thickest part of the slope body before the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine is 20 m, and the fully weathered layer and strongly weathered layer of Jurassic Yan'an group mudstone are 10-20 m thick.the Tanshan Coal Mine landslide is a bedrock landslide, and the main body of the landslide is Jurassic Yan'an group mudstone.The influence of mudstone on slope stability mainly includes three aspects: disintegration, expansibility and low strength.Disintegration refers to the nature of disintegrating dispersion of rock and soil mass after immersion in water, which is due to the unbalanced growth rate of the diffusion layer due to the uneven entry of water into voids and cracks, resulting in local repulsion greater than suction, and rock and soil mass collapse along the disintegration surface.Mudstone belongs to clay rock and has fine particles and poor engineering properties.It easily disintegrates and cracks into fragments or flakes upon contact with water.For example, colluvium is common at the slope toe of exposed mudstone free face; the drill core slowly cracks into fragments over time.In addition, because the mudstone contains hydrophilic clay minerals such as montmorillonite and illite, the mudstone expands in volume after water absorption and generally has expansibility, and the strength of the expansibility increases with the increase in montmorillonite content.At the same time, due to disintegration of the mudstone, the disintegration surface increases the specific surface area of mudstone and increases the water absorption capacity and expansibility of mudstone.The low strength of mudstone is closely related to its mineral composition, rock structure and water content.The mineral composition of mudstone is mainly soft hydrophilic minerals such as montmorillonite and illite, which have argillaceous cementation with naturally low compressive strength; in terms of rock structure, it is mainly affected by particle size, porosity and void distribution.Generally, rocks with large porosity, that is, rocks with poor compactness, have lower compressive strength; the influence of water content on mudstone strength is also very significant, and strength generally decreases with increased water content.The Jurassic Yan'an group mudstone contained in the pre-sliding slope of the Tanshan Coal Mine contains more clay minerals, with strong hydrophilicity, prominent hydration, significant disintegration and expansion, low strength and poor engineering properties, which provide a material basis for the formation of the Tanshan Coal Mine landslide [54].

Inducing Factors (1) Coal mining activities
The pre-sliding slope gradient of the Tanshan Coal Mine was only 11 • .Through investigation of geological disasters in the Yuanzhou District, Guyuan City, Ningxia, sliding failure generally does not occur on such a small slope, while the formation of the Tanshan Coal Mine landslide is closely related to human engineering activities.The main engineering activities in the area include slope cutting, house building, cultivated land reclamation, highway plant construction, coal mining, etc., and the most influential factor for the formation of the Tanshan Coal Mine landslide is the influence of the Tanshan Coal Mining Subsidence Area.
The deformation and failure of the overlying rock and soil mass are caused by underground mining.The coal mining subsidence area formed by underground coal mining changes the original stress field of the slope body, resulting in the destruction of the surrounding rock and soil mass, the outward transmission of displacement and the formation of collapse basin and cracks on the surface.Due to the action of joints and fissures in rock and soil, its movement and deformation is discontinuous and asymmetric.The front slope is displaced towards the coal mining subsidence area under the influence of the coal mining subsidence area, which leads to the loss of support for the rear slope and creep deformation in the direction of the free face under the action of self-weight stress [11,55].The movement and deformation of rock stratum destroys the stable state of the previous slope and reduces the mechanical strength of the weak structural plane in the rock mass so as to induce the landslide together with the infiltration of surface cracks and water.
(2) Rain The rainfall in the study area is characterized by a concentrated rainy season and uneven precipitation distribution.From July to September, precipitation accounts for more than half of the annual precipitation, and most is sudden regionally heavy precipitation, which is an important inducing factor for the frequent occurrence of geological disasters such as landslides in this area.Its impact on slope stability is mainly reflected in the following three aspects: 1 It reduces the strength of rock and soil [43].Increasing the water content of rock and soil concurrently reduces the effective cohesion c and effective internal friction angle ϕ , which directly leads to the reduction of the shear strength of the rock and soil mass. 2 It Increases the unit weight of rock and soil [56].When rainfall is sufficient, rainfall and surface water infiltrate joints and fissures, increasing the unit weight of the rock and soil mass.From the limit equilibrium method of slope stability calculation, it can be seen that the increase in rock and soil unit weight increases the sliding force on the one hand and the normal force on the other hand so as to increase the anti-sliding force.However, it can be calculated that the increase in the sliding force is always greater than that of the antisliding force, which eventually leads to the reduction of the safety factor. 3It Influences fracture water pressure and buoyancy pressure [57].Precipitation infiltrates along joints and fissures, increasing fissure water pressure and the further deepening and penetration of fissures.At the same time, precipitation infiltration may raise the groundwater level and increase the buoyancy force, which reduces the stability of the slope.In addition to the above three aspects, the impact of surface water on slope stability also plays an important role in scouring the slope toe and slope surface, changing the geometric shape of the slope, forming free faces, and providing spatial bases for the formation of a landslide.
In the formation process of the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine, precipitation and surface water are the main reasons for the new deformation cracks and resurrection trend of the landslide.The gully at the front edge of the slope scours the slope toe to form a steep free face, which provides space for the slope to slide towards the slope toe.At the same time, precipitation and surface water infiltrate along the fissures formed in the coal mining subsidence area, forming a weak zone at the contact surface of strongly and moderately weathered mudstone, reducing the shear strength of the rock mass and finally leading to the formation of a landslide in the Tanshan Coal Mine.

Analysis of Formation Mechanism of Overlying Landslide in the Coal Mining Subsidence Area
The mining of underground minerals may induce secondary geological disasters such as collapse, landslide, debris flow and ecological environment deterioration, resulting in heavy casualties and economic losses.Therefore, the study of the formation mechanism of landslides in coal mining subsidence areas is of great significance to predict and prevent similar landslides [33].According to previous numerical simulation results, under the influence of only the coal mining subsidence area, the collapse basin and mining fractures are far from the critical state that can cause slope slide and are inconsistent with the field investigation.It can be judged that rainfall also has a great impact on the landslide stability of the Tanshan Coal Mine.The Tanshan Coal Mine landslide was formed in the rainy season of 1966.Underground mining destroyed the stress field of the slope, resulting in tension cracks at the rear edge of the slope and collapse basins and mining cracks at the front edge.Rainfall and surface water infiltrated downward along the rock cracks formed by tension, mining and coal mining subsidence area collapse.Due to the low strength and high permeability of strongly weathered mudstone compared with moderately weathered mudstone, after rainfall and surface water infiltration, it is easy for water to combine with mudstone at the contact surface between moderately weathered mudstone and strongly weathered mudstone, resulting in mudstone water saturation and softening, reducing the shear strength of rock and soil mass and forming local weak zones; at the same time, due to the continuous occurrence of coal mining subsidence area collapse, the cracks in the rock and soil mass further increase, and the process continues to go back and forth until local softening finally forms a weak zone through the material.At the same time, the self-weight of the upper rock and soil mass increases after water absorption.Under the joint action of these factors, the sliding failure of the slope finally occurs.The whole sliding process of the Tanshan Coal Mine landslide can be divided into four stages [8,33].

Slope Deformation Stage
At this stage, there is a small collapse failure in the coal mining subsidence area, the rock and soil stress field causes movement of the rock stratum, the stress state of the rock mass changes, the rock stratum above both ends of the coal mining subsidence area is compressed, and the support pressure is transferred to the depths on both sides; the roof and floor of the working face in the coal mining subsidence area form a decompression zone, and the rock stratum in the lower part has undergone small elastic recovery deformation.Under the action of gravity, the roof rock layer first bends towards the coal mining subsidence area; then it breaks, falls and sinks in layers; and then separation, cracks and fractures occur.

Slope Creep Stage
At this stage, the coal mining subsidence area completely collapses, the direct roof completely loses its strength, bends to the coal mining subsidence area and destroys it, and then breaks into rock blocks of different sizes to fall down and fill the coal mining subsidence area.The middle zone strata in the middle and upper part of the full coal mining subsidence area will sink in layers under the action of self-weight stress and the lower strata, and the falling broken rock blocks in the coal mining subsidence area are gradually compacted.The displacement of the rock stratum is gradually transmitted to the surface, forming a collapsed basin on the surface.Due to the non-uniformity of slope movement and deformation, mining cracks are generated in the surface tensile deformation area outside the coal mining subsidence area.Mining cracks may not only crack the rock mass and directly affect the strength of the rock mass, but they also aggravate the infiltration of rainfall and surface water, further reducing the strength of the rock mass and affecting the stability of the mining mountain.
The mining cracks caused by the coal mining subsidence area and the tensile cracks formed at the trailing edge of the slope's movement accelerate the infiltration of rainwater, which combines with mudstone at the contact surface of strong weathering and moderate weathering to reduce the strength of mudstone and form a weak zone, and the slope produces creep deformation along the weak zone.In this process, the upper rock mass slides under the action of self-weight stress, the tension crack at the trailing edge is further widened and deepened, the sliding body is gradually separated from the trailing edge, the vertical displacement and horizontal displacement increase, and a small-scale collapse has occurred on the free face of the leading edge.
According to the survey, before the landslide was formed in 1966, local residents found tension cracks at the rear edge of the slope, accompanied by small beaded sinkholes.When it rains, rainwater seeps along the cracks and sinkholes.There are signs of collapse at the front and mining fractures at the edge of the collapsed basin.This shows that the stress field of the slope has changed due to the influence of the coal mining subsidence area.The increase in rock and soil cracks and the expansion of tension cracks are conducive to rainfall and surface water infiltration, which provides a prerequisite for the next sliding mass movement.

Landslide Movement Stage
After a long-term creep deformation stage, the gob collapse gradually reduces the supporting effect of the lower rock mass on the upper rock mass, the tensile cracks at the rear edge of the slope continue to increase, and the rock mass strength continues to decrease.At the same time, combined with the infiltration of rainwater and surface water, the strength of the weak zone decreases and runs through into a continuous sliding surface.At the same time, the mining cracks in the front also increase, which not only reduces the strength of the lower rock mass, but also increases the channel of rainwater infiltration and increases the gravity of the rock mass.When the sliding force exceeds the anti-sliding force provided by the sliding bed, the sliding body obtains a relatively large initial velocity under the action of self-weight stress, and the potential energy of the sliding body is continuously transformed into kinetic energy, so the sliding body slides violently along the sliding surface and cuts out along the shear outlet and slides to the front gully to form a landslide.

Landslide Accumulation Stage
After violent sliding and energy release, the sliding speed of the Tanshan Coal Mine landslide decreases slowly, and finally the sliding mass accumulates at the front edge of the landslide and the gully in front of the slope.

Comparison of Calculation Methods of Landslide Stability Coefficient
At the same time, it can be seen from Table 5 that there are obvious differences between the calculation results of the limit equilibrium method and the finite element method.Take the M-P method and finite element method in the limit equilibrium method under the rainfall state as an example.
The difference in the stability coefficient calculated by the two methods is about 5%, which is mainly related to the normal stress distribution along the sliding surface.Stress concentration occurs at the slope toe of the landslide mass, and the normal stress at the bottom of the strip in the limit equilibrium method is mainly caused by the gravity of the strip.Therefore, the limit equilibrium method cannot obtain local shear stress concentration, which is a limitation of the limit equilibrium method.Comparatively, the stress calculated by the finite element method is closer to the real situation, and one advantage of the finite element method is that the local safety factor of each soil strip can be observed.Figure 14 shows the change in the stability coefficient from the first strip to the last strip of the sliding surface obtained by the finite element method.At the same time, it can be seen from Table 5 that there are obvious differences between the calculation results of the limit equilibrium method and the finite element method.Take the M-P method and finite element method in the limit equilibrium method under the rainfall state as an example.
The difference in the stability coefficient calculated by the two methods is about 5%, which is mainly related to the normal stress distribution along the sliding surface.Stress concentration occurs at the slope toe of the landslide mass, and the normal stress at the bottom of the strip in the limit equilibrium method is mainly caused by the gravity of the strip.Therefore, the limit equilibrium method cannot obtain local shear stress concentration, which is a limitation of the limit equilibrium method.Comparatively, the stress calculated by the finite element method is closer to the real situation, and one advantage of the finite element method is that the local safety factor of each soil strip can be observed.Figure 14 shows the change in the stability coefficient from the first strip to the last strip of the sliding surface obtained by the finite element method.It can be seen that the stability coefficients of each strip obtained by the finite element method are different, and the difference between the maximum value and the minimum value can reach 4.6.Combining the available anti-sliding shear forces and sliding shear forces on the whole sliding surface and averaging them, a stability coefficient similar to that provided by the limit equilibrium method can be obtained.It can be seen from Figure 14 that the stability coefficient of the strip is closely related to the shape of the sliding surface.The trailing edge of the sliding body possesses a high safety factor due to the small division of the blocks, the small self-weight of the blocks, and the lesser sliding force relative to the shear strength; in the middle of the landslide, the stability coefficient decreases obviously at the junction of the landslide platform; at the landslide tongue at the front edge of the sliding mass, the angle of the sliding zone is very small and rock and soil mass accumulate here, so the safety factor increases obviously; finally, at the edge of the gully, the safety factor is very low due to the sharp increase in the angle of the sliding band and the influence of the free face.
From the above analysis, it can be seen that the use of finite element calculation in It can be seen that the stability coefficients of each strip obtained by the finite element method are different, and the difference between the maximum value and the minimum value can reach 4.6.Combining the available anti-sliding shear forces and sliding shear forces on the whole sliding surface and averaging them, a stability coefficient similar to that provided by the limit equilibrium method can be obtained.It can be seen from Figure 14 that the stability coefficient of the strip is closely related to the shape of the sliding surface.The trailing edge of the sliding body possesses a high safety factor due to the small division of the blocks, the small self-weight of the blocks, and the lesser sliding force relative to the shear strength; in the middle of the landslide, the stability coefficient decreases obviously at the junction of the landslide platform; at the landslide tongue at the front edge of the sliding mass, the angle of the sliding zone is very small and rock and soil mass accumulate here, so the safety factor increases obviously; finally, at the edge of the gully, the safety factor is very low due to the sharp increase in the angle of the sliding band and the influence of the free face.
From the above analysis, it can be seen that the use of finite element calculation in the framework of the limit equilibrium method has obvious advantages.It not only overcomes the numerical convergence difficulties of the traditional limit equilibrium method in some cases, but it is also more in line with the actual stress distribution, making the calculation results more accurate and reliable.

Conclusions
Based on the detailed investigation data of geological disasters in the Yuanzhou District, we paper systematically study the formation mechanism and stability of the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine by means of geotechnical engineering exploration, indoor geotechnical tests, theoretical analysis and discrete element finite element numerical simulation, and we draw the following main conclusions: (1) The slightly inclined bedding mudstone landslide in the Tanshan Coal Mine is 850 m long from east to west, 500 m wide from north to south, 10-20 m thick, 475,000 m 2 in total area, 10,875,000 m 3 in volume and has a sliding direction of 295 • .The main body of the landslide is composed of Jurassic mudstone.From the material composition and scale of the landslide mass, the landslide belongs to a giant bedrock landslide.From analysis of landslide deformation characteristics and movement properties, the landslide is a traction landslide caused by a coal mining subsidence area.(2) Discrete element simulation analysis considering only the coal mining subsidence area shows that when the coal mining subsidence area collapse is completed, the creep deformation of the slope is not enough to lead to the formation of a landslide, which indicates that the formation of the landslide in the Tanshan Coal Mine is closely related to rainfall and surface water infiltration.When the contact surface between strong weathering and moderate weathering is set as a weak zone in the model, the simulation results are consistent with the actual sliding situation of the landslide.(3) The formation of the slightly inclined bedding mudstone landslide in the Tanshan Coal Mine is mainly due to the joint action of internal factors and inducing factors.The internal factors are mainly controlled by geotechnical types and engineering geological properties.The inducing factors are mainly coal mining activities and rainfall.(4) The sliding process of the landslide is divided into four stages: slope creep, slope deformation, landslide movement and landslide accumulation.(5) By comparing the calculation results of the limit equilibrium method and the finite element limit equilibrium method, we found that under the rainfall condition, the stability coefficient obtained by the two methods is about 5%, which is mainly due to the fact that the normal stress at the bottom of the strip in the limit equilibrium method is caused by the gravity of the strip and cannot simulate the actual stress distribution.The finite element limit equilibrium method not only calculates the stress closer to that of the real situation, but it also provides the local safety factor of each soil strip.Therefore, selecting the appropriate constitutive model, using the finite element method to simulate the stress distribution of the slope, and then importing it into the limit equilibrium equation for calculation makes the calculation results more accurate.

Figure 1 .
Figure 1.The shape of landslide in the Tanshan Coal Mine.

Figure 2 .
Figure 2. Geological map and section of landslide in coal mining subsidence area of Tanshan, China.

Figure 1 .Figure 1 .
Figure 1.The shape of landslide in the Tanshan Coal Mine.

Figure 2 .Figure 2 .
Figure 2. Geological map and section of landslide in coal mining subsidence area of Tanshan, C Figure 2. Geological map and section of landslide in coal mining subsidence area of Tanshan, China.

Figure 3 .
Figure 3. Discrete element model of the slope before the landslide in the Tanshan Coal Mine.

Figure 3 .
Figure 3. Discrete element model of the slope before the landslide in the Tanshan Coal Mine.
lower boundary is X-and Y-direction displacement constraints, and the upper boundary is a free boundary.The model is taken as an approximate ideal elasticplastic model, and the Mohr-Coulomb yield criterion is selected as the failure criterion.

Figure 4 .
Figure 4. Finite element calculation model of slope before landslide and after mining in the Tanshan Coal Mine.

Figure 4 .
Figure 4. Finite element calculation model of slope before landslide and after mining in the Tanshan Coal Mine.

Figure 5 .
Figure 5. Numerical simulation results under the coal mining subsidence area.

Figure 5 .
Figure 5. Numerical simulation results under the coal mining subsidence area.

Figure 10 .
Figure 10.Displacement vector of slope before sliding.

Figure 10 .
Figure 10.Displacement vector of slope before sliding.

Figure 10 .
Figure 10.Displacement vector of slope before sliding.

Figure 11 .
Figure 11.Displacement contour map of the slope in the X-direction before sliding.

Figure 11 .
Figure 11.Displacement contour map of the slope in the X-direction before sliding.

Figure 12 .
Figure 12.Y-direction displacement contour map of slope before sliding.Figure 12. Y-direction displacement contour map of slope before sliding.

Figure 12 .
Figure 12.Y-direction displacement contour map of slope before sliding.Figure 12. Y-direction displacement contour map of slope before sliding.

Figure 13 .
Figure 13.The contour map of major principal stress after mining.

Figure 13 .
Figure 13.The contour map of major principal stress after mining.

Table 1 .
Calculation parameters of discrete element numerical simulation before landslide in the Tanshan Coal Mine.

Table 2 .
Calculation parameters of finite element numerical simulation before landslide in the Tanshan Coal Mine.

Table 3 .
Differences in limit equilibrium methods.

Table 4 .
Stability calculation parameters after landslide in the Tanshan Coal Mine.

Table 5 .
Calculation results of landslide stability coefficient after sliding in the Tanshan Coal Mine.