Discrete element analysis of deformation features of slope controlled by karst fissures under the mining effect: a case study of Pusa landslide, China

Abstract Karst landforms are widely distributed in the southwestern mountain areas of China, and the continuous underground mining activities lead to frequent occurrence of catastrophic collapses and landslides. Revealing the relationship between the development characteristics of the controlling karst fissures and the slope deformation process is crucial to understand the collapse and landslide phenomena. The Pusa landslide is selected as the geological prototype of discrete element analysis, and the universal distinct element code (UDEC) is applied to simulate the overall deformation response of the mountain containing extensive karst fissure during the mining process. The results show that under the action of mining, the roof above the goaf bends and subsides, and the middle of the roof even breaks and collapses. The separation fractures effectively block the upward transmission of the collapse state of the rock stratum. The bottom of the karst fissure is susceptible to cracking first in the process of coal seam mining due to stress concentration, and the area of severe deformation in the slope coincides with the mining pressurization area. The morphology of the karst fissure controls and determines the deformation characteristics of the rock mass at the slope top, and only the karst fissure located within the mining influence range is the object to be considered in the slope stability analysis. The limit karst fracture depth, about 1/3 of the slope height, is the limit value to determine whether the rock mass at the slope top is toppled or slipped. The relationship between the karst fissure and the free surface gradually changes from the directional or co-directional to the reverse, the motion state of the rock mass at the slope top changes from slipping to toppling, and the role of karst fissure changes from a potential slip surface to the cracking boundary. Although the deformation damage of the reverse structural slope is not very serious, the influence of the karst fissure on the stability of the slope still cannot be ignored. This study aims to provide basic theoretical support for the subsequent research on the failure mechanism of karst mountains under the combined action of multi-structural planes.


Introduction
The two-dimensional planar geological interface with certain direction, scale and form existing in the rock mass under tectonic stress, including material differentiation surface and discontinuity surface, such as bedding plane, foliation, faults, joints, cleavages, fissures and other fracture systems, is collectively referred to as structural plane (Peacock 2001;Richard et al. 2008;Peacock et al. 2018).The presence and abundance of discontinuities in the rock mass are known geological factors controlling the stability of slopes (Gu 1979;Miller 1983;Guzzetti et al. 1996;G€ unther 2003;Fourniadis et al. 2007;Marchesini et al. 2015), and determine the formation and volume of rock falls, topples, rock slides and rock avalanches (Cruden 2003;Katz and Aharonov 2006;Santangelo et al. 2015;Henriques et al. 2015;Esposito et al. 2021).
Generally, large rock landslides are divided into two categories: structural planecontrolled landslides and non-structural plane-controlled landslides (Hungr and Evans 2004).The structural plane cuts the rock mass into blocks of different macroscopic scales, seriously destroying the integrity of the rock mass (Zhao et al. 2019;Guo et al. 2020).More importantly, the structural plane is often combined with weak rock strata and slope surface to form potential slip surface or tension crack boundary, thus contributing to the occurrence of large-scale landslides or collapses (Swanson 1992;Schultz 2000).The collapse-debris flow or landslide-debris flow disasters induced by underground mining are characterized by fast movement, wide disastercausing range, huge destructive force and difficult prediction (Yang et al. 2015;Ge et al. 2018;He et al. 2019).The process of deformation and failure of mountain containing extensive structural planes is often accompanied by loosening, extension and cracking changes of structural planes.Underground mining even changes the original features of the structural plane, especially the controlled weak structural plane (Li et al. 2016;Ma et al. 2018;Fern andez et al. 2020).
The large-scale collapse and slide disasters in karst mountainous areas of the southwest China mostly occur in the wing or core of the folded structural belt composed of layered carbonate rocks and clastic rocks.The dissolution of carbonate rocks often leads to a large number of karst fissures and cavities in this kind of hard rock mass, and the rock strength of key blocks deteriorates severely under the effect of seasonal seepage of karst fissures (Liu et al. 2020a).In landscapes dominated by layered rocks, defining the geometrical relationships between the attitude (strike and dip) of the controlling structural plane and the attitude (gradient and aspect) of the terrain is crucial to understand landslide phenomena (Qi et al. 2004;Guzzetti et al. 2008;Santangelo et al. 2015).Song (2022) argues that in the study of slope structure, the influence of the combination relationship between the controlling structural plane and the weak rock strata or slope surface on the slope stability should be considered, which is of enlightening significance for relevant research.The results of current studies on the influence of controlling structural plane on slope stability are mostly qualitative evaluation, while the quantitative analysis on the relationship between two-dimensional morphological characteristics of structural plane and slope deformation features is less.Moreover, the relationship between underground coal mining and slope deformation process controlled by structural plane has not been clearly explained.
This study first summarizes the basic characteristics of typical mining-induced landslides and collapses in karst mountain area in the southwest China since 1980, and emphasizes the importance of discontinuities inside the mountain in controlling the occurrence of disasters.Then, the case of '8.28' mega landslide in Pusa, Nayong County, Guizhou Province is selected as the research object, and the mountain with large deformation degree at the late stage of mining (from 2013 to 2017 before the collapse accident occurred) is generalized into a two-dimensional model, which contained karst fissure of different lengths (development depth), dips and position (distance from the slope surface).The discrete element simulation method is applied to reproduce the complete process of deformation and fracture evolution of the mountain containing extensive karst fissures during the mining process.By analyzing the strata movement characteristics and slope stress distribution law under the action of mining, the relationship between karst fissure morphology and slope deformation features is revealed.Considering that multiple karst fissures may have superimposed effects on the deformation and failure of slopes, this paper only discusses with a single karst fissure to get the corresponding basic conclusions.

Study area
The Yunnan-Guizhou Plateau and the surrounding areas in the southwest China are the largest contiguous exposed karst area in the world, covering an area of 6.2 Â 10 5 km 2 , accounting for about 1/15 of the national territory.The region is the base for China's high reserves of high-grade mineral resources and abundant hydro energy resources, as well as an important output site for west-east power transmission (Liu et al. 2015;Chen et al. 2021a).The southwest karst mountainous area is located in the Upper Yangtze platform and has formed a unique strongly folded geomorphology through multiple phases of tectonic movements, shaping steep and broken geomorphic landscape of fenglin karst and fengcong karst (Marjorie 1995;Waltham 2008).The surface of the mountainous area is covered by thick carbonate strata of the Permian and Triassic System, forming a complex regional geological environment with a strata structure by hard tops and soft bottoms or soft-hard interbedding, topography by its steep tops and gentle bottoms (Xu et al. 2016;Chen et al. 2021b).The large surface incision and topographic slope make the mountain body present the morphological characteristics of nearly upright in the upper part and gentle in the lower part, which provides huge potential energy for the formation of high-speed and long-distance debris flow after the instability of the mountain (Fomenko et al. 2021).According to the lithology classification standard of geology, the strength index is the basis for the classification of soft and hard rock strata.A class of rocks with uniaxial compressive strength (r c ) between 0.5 Mpa and 25 Mpa is defined as soft rock by International Society of Rock Mechanics (ISRM) (Hudson and Harrison 1997).Soft rock mainly includes shale, mudstone, sandstone, siltstone and so on, which has the characteristics of low strength, large porosity, low cementation and obvious weathering.In China, 29.4 MPa is commonly used as the limit value of uniaxial compressive strength for soft rock and hard rock (Yang 1986).Hard rock mainly includes quartzite, marble, granite, limestone and so on.
Subject to long-term dissolution, complex karst pipeline system had been formed inside the mountain, which provided advantageous channels for rainwater infiltration and induced concentrated seepage in the lithological contact zone.A series of karst fissures, karst pipes, karst caves and undulating folded geomorphological, resulting in the extremely uneven spatial and temporal distribution of groundwater, as well as frequent and substantial changes in hydrodynamic conditions, significantly deteriorate the stability of the slope (Guti errez et al. 2014).Unlike the mining of thick or medium-thick ore layers under the loose clayey soil layers and soft sediments in plain areas, the southwest mountains are mainly mined for thin ore layers, which are covered by thin layers of loose clayey soil, weak rock layers and thick layers of hard rock (Li et al. 2020a).The mining activities of underground minerals, while satisfying human energy needs, also change the initial stress environment of the rock strata in the goaf and the engineering geological conditions of the mining area (Lana 2014;Salmi et al. 2017).It not only causes surface subsidence and building deformation, but even triggers landslide disasters with volumes ranging from tens of thousands to tens of millions of cubic meters, resulting in serious human and property losses (Liu et al. 2019;Akcin 2021).Table 1 statistics the typical mining-induced landslides and collapses in karst mountain areas in the southwest China since 1980, focusing on the geological conditions, slope structure, influencing factors and deformation mechanisms for the occurrence of the disasters.Our research team has made detailed field investigation and research on the disaster characteristics, geological background and failure mechanism of Jiguanling landslide (Li et al. 2016), Zengziyan collapse (He et al. 2019), Jiweishan landslide (Li et al. 2021) and Pusa collapse (Yang et al. 2020;Yang et al. 2022a).Relevant information on the remaining geological disasters is obtained through literature collection.It is found that long-term and large-scale underground mining activities accelerate the deformation and stress adjustment of the mountain, resulting in tensile-shear fracture and compression-shear slip of the hard structural plane, and eventually cause the deformation and instability of the mountain (Chen et al. 2020;Wang et al. 2020a;L opez-Vinielles et al. 2020;Yang et al. 2022b).
Pusa collapse is a typical remote chain avalanche disaster induced by underground mining, which occurred in the back hill of Pusa village, Zhangjiawan town, Nayong county, Guizhou province, China (N26 38 0 04.55 00 , E105 26 0 56.14 00 ) (Figure 1).A small amount of cracks and small collapses occurred in the mountain during 2009-2015, and the mountain began to deform at an accelerated rate in 2016.On August 28, 2017, a rock mass of with the volume of 49.1 Â 10 4 m 3 descended at a high level, scraping off and accumulating loose deposits to finally form an 82.3 Â 10 4 m 3 mass that destroyed the residential area at the slope foot, resulting in the death of 26

Model establishment and parameter determination
Taking the Pusa collapse as the geological prototype for numerical simulation research, a two-dimensional discrete element numerical model of the mountain before collapse is established along the main sliding direction at a scale of 1:1, as shown in Figure 3.In accordance with the test and calculation requirements of the technical specification 'standard for test methods of engineering rock mass' (GB/T 50266-2013), laboratory tests are carried out on various types of rock and coal seam samples collected on site.The density, bulk modulus, shear modulus, cohesion, internal friction angle and tensile strength of various materials are obtained by density test, uniaxial compression test, triaxial compression test and uniaxial tensile test, collated in Table 2.The determination of parameters for joints and rock interfaces is based on the field survey and engineering analogy method (Table 3), the main references are Zhong et al. (2020), Xiong et al. (2021), andCui et al. (2022).Universal Distinct Element Code (UDEC) is a two-dimensional numerical program based on the distinct-element method for discontinuum modeling.UDEC simulates the deformation response of a discontinuous medium, represented as a combination of discrete blocks, subjected to either static or dynamic loading.It has been widely used in the study of progressive failure of rock slope, permeability characteristics of rock mass, and sizetime effects of rock strength (Israelsson 1996).A variety of post-processing graphics can be displayed and output in UDEC, such as stress nephograms, displacement nephograms, joint state nephograms, seepage nephograms.It is also possible to set up different output items and extract specific values for each data in a certain range (points, lines and surfaces) before the nephogram is output.Coal and rock materials are regarded as elastoplastic materials and Mohr-coulomb plasticity failure criteria is adopted; the joints and the discontinuity surface between the rock strata are applied with the point contact-coulomb slip model in the UDEC.During the numerical calculation, the gravity g ¼ 9.8 N/kg is set.
The design of the two-dimensional numerical model has gone through two stages: preliminary conceptual design and logical design with necessary simplifications based on the test requirements, corresponding to Figure 3 and Figure 4, respectively.The bottom boundary of the numerical model in Figure 3 is 600 m long.It can be    In addition, it is found by collating and summarizing (Table 1) that the original slope line inclination of karst mountain source area in the southwest China is mostly concentrated in the range of 45 to 90 , while the original slope line inclination of the deposition area is from 15 to 40 .In order to more intuitively compare the deformation characteristics and fracture evolution process of slopes under different working conditions, the numerical model of Figure 3 is further simplified.It is important to note that the necessary simplification of the numerical model must follow certain design logic as well as reasonable assumptions.The free surface boundary of the model is all set to a straight line, the top and bottom of the slope are 325 m and 123 m high, respectively.The gradient is 70 , and the length of the slope projection to the horizontal plane is 96 m (from the horizontal coordinate 295 m to 391 m) (Figure 4).Three sets of structural planes are mainly developed in the rock mass of the collapse source area, with the spacing of 3 $ 8 m.In order to facilitate modeling, the joints occurrence of the two-dimensional model are 310 /7 , with joints spaced at 4 m to 6.2 m and the orthogonal secondary joints at 8 m to 12 m.Only one karst fissure is set at the shoulder of the slope with a maximum width of 2 m, and no fault is set below the steep slope of the source area.

Layout of the monitoring points
In order to reveal the characteristics and mechanisms of the dynamic response of karst mountains containing deep and large fissures under the action of mining, typical monitoring lines and data points are set up inside the slope model for analysis (Figure 5).Four lateral displacement monitoring lines, numbered A, B, C, and D, are selected from top to bottom.Lines A and B traverse the karst fissures, and lines C and D are located below the karst fissures.A speckle is selected as the data point at a certain interval on the monitoring line, and a total of thirteen columns of the data points are divided.Besides, a monitoring line (E) is also set along the slope surface, with ten data points set at equal intervals from top to bottom.By counting the horizontal or vertical displacement of all data points, the displacement variation curves of the five monitoring lines during the mining process are obtained.The analysis of displacement variation is shown in Sections 4.1.2,4.2.2 and 4.3.2.

Coal seam excavation scheme
The mining of Pusa coal mine adopts the inclined well development and strike long wall coal mining method, and the roof management of goaf adopts the fallen method.The layout direction of the working face is consistent with the inclined direction of the coal seam.This mining method is easy to cause the roof of the coal seam to collapse due to large area exposure.The storage situation of coal seam and mining range of Pusa coal mine are shown in Figure 2. The main coal seams in the study area in the early stage of mining (from 1995 to 2012) were M20 and M16 (Figure 2).The deep distribution of coal seams, the wide range of mining activities and the uneven distribution of excavation areas led to a low level of mountain deformation, so the mountain was difficult to be generalized to a two-dimensional model.At the late stage of mining (from 2013 to 2017 before the collapse accident occurred), the main coal seams in the study area were M14 and M10, with M6 presumably being partially mined (Figure 2).Zhong et al. (2020) believe that there is a certain time correspondence and potential causal relationship between the tension cracks behind the source area and the mining of M14 and M10.Li et al. (2020b) highlights that the breakage of thick, hard intermediate strata of multiple coal seams has a decisive influence on the movement of overburden, while soft, thin intermediate strata are deformed simultaneously with the overburden.The thicknesses of M6, M10, M14, M16, and M20 are 2.01 m, 2.12 m, 1.23 m, 1.49 m, and 0.84 m, respectively, and the total thickness of the coal seam is 7.69 m.The intermediate strata in M14 and M10 are silty mudstone, around 15 m thick, in line with the characteristics of soft, thin intermediate strata.A single coal seam (M14) is thus used instead of multiple coal seams in the mountain model, with the thickness of 8 m, and the effect of repeated mining was not considered.The length of mining working face is 200 m, and no coal pillar is left during excavation.The single excavation step is 20 m, and the excavation is carried out sequentially along the mining direction, with a total of ten excavations (Figure 4a).The model is in the initial stress equilibrium before coal seam mining, and the equilibrium calculation of the model is performed after each excavation.When the model is in the stress equilibrium state, the next excavation is carried out.When the maximum unbalanced force in the model tends to infinity or the ratio of the current maximum unbalanced force to the initial maximum unbalanced force is less than 10 À5 , or the nodal displacement tends to a constant and the value no longer changes, the model can be determined to be in a stress equilibrium state.

Simulation of operating conditions
Extensive karst fissures are a kind of two-dimensional planar geological interfaces with directionality, and the cross-section is wedge-shape, extending gradually from the top of the slope to the interior of the slope.As can be seen from the section of Pusa landslide along the main sliding direction (Figure 2), the two-dimensional morphological characteristics of the karst fissures are mainly described by three parameters: length (development depth), dip and position (distance from the slope surface).According to the actual distribution characteristics of karst fissures in the mountain, the development depth of karst fissure, the distance from the bottom of the karst fissure to the slope surface, and the dip of karst fissure are all taken as the research variables.The depths of karst fissure are set as 0 m, 40 m, 70 m, and 100 m, and the distances from the bottom of the karst fissure to the slope surface are 60 m, 90 m, and 120 m, respectively.Since karst fissures are mostly developed in the near-vertical directions, four dips of karst fissures are set in this study, which are 70 , 80 , 90 , and 100 , respectively (Figure 4).In this paper, a total of eight operation conditions (Table 4) are set up, and  the two-dimensional calculation model for each operation condition is shown in Figure 4, where operation condition 0# is a blank control model.

Discrete element simulation results
4.1.Relationship between karst fissure depth and slope deformation features deformed area and un-deformed area, and the deformed area can be further subdivided into three parts: the roof above the goaf, the rock mass near the slope surface, and the rock mass at the slope top (Figures 6a-d and 7a-d).The middle of the roof may break and collapse due to excessive bending deflection (Liao et al. 2020).Except for the horizontal movement of the rock strata at the left fracture boundary of the roof towards the inside of the slope, both the middle of the roof and the rock strata at the right fracture boundary move horizontally towards the free surface (Figure 7a-d).The emergence of separated fractures effectively blocks the upward transmission of rock strata collapse state (Wu et al. 2018;Ning et al. 2020), making the vertical deformation of the rock mass above the separated fractures decrease significantly, while the horizontal deformation gradually decreases from the free surface toward the inside of the slope.The settlement centre line of the slope does not coincide with the midpoint line of the goaf, but shifts to the left side, which is related to the inclination of the coal seam (Figure 6a-d).In general, the vertical and horizontal deformation features of the roof and the rock mass near the slope surface are basically the same in different slopes containing karst fissure of different depths; however, the deformation features of rock mass on the left side of the karst fissure show obvious differences.
As the development depth of the karst fissure increases from 0 m to 70 m, the vertical displacement of the rock mass on the left side of the karst fissure gradually increases, and the vertical displacement from the karst fissure to the slope surface gradually increases within the region (Figure 6a-c).Similarly, the horizontal displacement of this part also increases with the increase of the depth of the karst fissure, while the trend of horizontal displacement from top to bottom gradually decreases within the region (Figure 7a-c).Through the above description of the displacement changes, it can thus be inferred that when the depth of the karst fissure is 0 m, 40 m and 70 m, the rock mass on the left side of the karst fissure tilts towards the free surface.Although the vertical displacement of the rock mass on the left side of the karst fissure increases when the depth of the structural plane is 100 m, however, the vertical displacement from the karst fissure to the slope surface gradually decreases within the region (Figure 6d).The horizontal displacement of this part decreases significantly, and the horizontal displacement increases slowly from top to bottom within the region.The horizontal deformation at the slope top is far less than that at the bottom of the karst fissure, indicating that the rock strata at the bottom of the karst fissure are extruded (Figure 7d).It can be speculated that when the depth of the karst fissure is 100 m, the rock mass on the left side of the karst fissure slides along the long karst fissure, resulting in the extrusion of the rock mass at the bottom of the karst fissure under the shearing effect of the sliding force.
From the above analysis, the development depth of the karst fissure can be used as an index to judge what kind of deformation process the rock mass on the left side of the karst fissure will undergo.By carefully observing the horizontal displacement nephogram of the slope without karst fissure (Figure 7a), it can be found that the deformation area near the slope surface is approximately triangular; while the deformation area near the slope top is approximately rectangular, and the horizontal displacement of the rock strata gradually decreases from top to bottom.Through measurement, the height of the rectangular deformation zone is about 65 m, about 1/ 3 of the slope height, and this distance can be considered as the limit karst fissure depth of the rock mass at the slope top for toppling failure.Once the depth of the karst fissure exceeds the limit depth, the karst fissure will act as a potential slip surface forcing the rock mass at the slope top to slide.
The vertical displacement of all data points on lines A to D and the horizontal displacement of all data points on line E are counted to obtain the displacement variation laws of five monitoring lines, as shown in Figure 8.The meanings of the negative values of displacement in the figures are the same as those in the displacement nephograms.
Under the conditions of different depths of karst fissure, the variation trends of vertical displacement in lines A and B are very similar.When there is no karst fissure in the slope, the vertical displacement of data points A1-A10 and B1-B10 gradually decreases from À0.25 m to À0.05 m, and the displacement curve is linearly distributed, indicating that the rock mass in the deep part of the slope is relatively more stable (Figure 8a and b).However, when there is a karst fissure with a certain depth in the slope, the vertical displacement of data points on the left side of the karst fissure decreases linearly; while the vertical displacements of the remaining data points on the right side of the karst fissure are basically maintained at À0.01 m (Figure 8a and  b).It can be seen that the karst fissure divided the slope into two discontinuous regions, which served to control the deformation development and block the failure transmission.Referring to the specific location of the monitoring line (Figure 5a), and the division criteria of the slope deformation area affected by mining, it is considered that lines C and D, located below the karst fissure, are in the area basically not affected by the karst fissure.The vertical displacements of both monitoring lines increase first and then decrease (Figure 8c and d), which is consistent with the trapezoidal or conical distribution characteristics of the displacement equipotential lines in the vertical displacement nephograms (Figure 6a-d).The vertical displacements of each data points on line D are most similar under the four conditions (condition 0# to 3#), indicating that the closer to the mining area, the less the deformation of the rock mass is affected by the karst fissure.When there is no karst fissure in the slope, the horizontal displacements of data points E1-E4 located near the slope top decrease gradually, and the horizontal displacements of the data points E5-E10 increase first and then decrease slightly.When the depths of the karst fissure are 40 m and 70 m, respectively, the horizontal displacements of the data points E1-E5 and E1-E6 corresponding to the karst fissure decrease linearly; subsequently, the horizontal displacements of E6-E10 and E7-E10 continue to decrease, but the change rate is relatively slow (Figure 8e).However, When the depth of karst fissure is 100 m, the horizontal displacement of line E increases first and then decreases (Figure 8e).Summarizing the horizontal displacement variation of line E can not only further corroborate the previous analysis of the motion state of rock mass on the left side of the karst fissure, but also reflect the deformation features of rock mass near the slope surface.Affected by rock strata inclination and coal seam mining, in addition to the downward bending deformation, the overburden of the goaf also undergoes lateral deformation due to the extrusion of the rock mass above (Wang et al. 2022).

Variation law of slope stress
The vertical stress nephograms of the initial equilibrium state of the slopes and the equilibrium state of the slopes after coal seam mining in conditions 0#, 1#, 2#, 3#, 4#, and 5# are shown in Figures 9 and 10, respectively.The negative value indicates compressive stress and the positive value indicates tensile stress.This topic referenced Figures 9a-d and 10a-d.
The initial vertical stress of the slope is layered distribution, gradually increasing from top to bottom, and the maximum compressive stress at the bottom of the model is about 7.80 MPa (Figure 9a).In the slope containing the karst fissure, the stress equipotential lines bulge upward on the left side of the bottom of the karst fissure, indicating that the compressive stress concentration occurred at the bottom of the karst fissure under the action of self-weight, and this area is also a fragile area prone to deformation and failure (Figure 9b-d).As the depth of the karst fissure increases from 40 m to 70 m and 100 m, the maximum compressive stress at the bottom of the karst fissure gradually increases from 2.75 MPa to 6.27 MPa and 11.70 MPa (Figure 9b-d), showing that the greater the depth of the karst fissure, the more obvious the concentration of compressive stress at the bottom of the karst fissure.
After the mining of working face, the roof of the goaf under all four conditions (condition 0# to 3#) has bending deformation, and the direct roof even breaks down.The interior of the slope tends to a new state of stress equilibrium.The left, middle and right sides of the deformed roof have formed stress-concentrated mining pressurization zone, among which the right side one is the largest in scope, extending obliquely downward from the bottom of the karst fissure to the coal seam in a striplike distribution (Figure 10a-d).The mining pressurization zone coincides with the area where the violent deformation occurred in the slope, and the karst fissure with sufficient depth becomes the dominant channel of slope deterioration, accelerating the process of slope destruction.

Fracture evolution process
The fracture evolution processes of the slope in conditions 0#, 1#, 2#, 3#, 4#, and 5# during coal seam mining are shown in Figure 11.Several grey blocks constitute the slope model in the figure, and the red line segments represent the open joint plane, referenced in this topic as Figure 11a-d.
During the mining process, the height, form and fissure distribution pattern of the deformation area above the goaf are basically the same in all conditions, but the fissure expansion at the bottom of the karst fissure or the slope top presents different characteristics with each forward excavation.After excavation of 20 m, the transverse cracks appeared in the direct roof above the goaf, and tension-shear cracks with oblique downward development appeared on the rock mass on the left side of the karst fissure only in condition 3# (Figure 11d -1).When a 60 m length was excavated, tension cracks appeared at the bottom of the karst fissure in condition 3#, and the already existing tension-shear cracks continued to develop downward (Figure 11d -2).As the working face advanced to 80 m and 100 m, tension cracks appeared at the bottom of the karst fissure in conditions 2# and 1# for the first time, respectively (Figure 11c -3 and b -4), while only the overburden above the goaf was still deformed in condition 0# (Figure 11a -4).According to the distribution characteristics of fissures, the deformation area of overburden above the goaf can be divided into three parts from bottom to top: caving zone, fractured zone and bending subsidence zone (Palchik 2015a, Palchik 2015b;Zhang et al. 2018).The caving zone was mainly distributed by separation fractures, and the middle part of the rock strata was prone to fracture.The fractured zone was dominated by the distribution of dense oblique stepped fissures, and the rock strata were subjected to shear failure.In the bending subsidence zone, short fissures were interconnected and the rock strata underwent bending deformation.The division of the above deformation area was in line with the division standard of 'three horizontal zones' of the overburden above the goaf proposed by Liu (1995).After excavation of 120 m, the tension crack appeared at 75 m from the slope top for the first time in condition 0# (Figure 11a -5), and the cracks at the bottom of the karst fissure continued to develop downward in other conditions (Figure 11b -5 to d -5)).Comparing the mining length of the coal seam corresponding to the first appearance of tension cracks at the slope top or the bottom of karst fissure in the four conditions, it can be found that the greater the depth of karst fissure, the earlier the appearance of tension cracks at the bottom of karst fissure (Figure 11a -5, b -4, c -3, and d -2).When excavated to 140 m, two extensive deformation areas in the slope were clearly discernible, which were the deformation area of the overburden above the goaf and the fracture extension area at the bottom of the karst fissure or the slope top (Figure 11a -6 to d -6).When a 200 m length was excavated, the roof above the goaf collapsed due to the bending deformation reaching the limit, and the mutual contact of two deformation areas increased the failure range of the slope (Figure 11a -7 to d -7; Egorov et al. 2001;Liu et al. 2020b).

Features of rock strata movement
To explore the relationship between karst fissure location and slope deformation features, it is necessary to ensure that the distance from the bottom of the karst fissure to the slope surface is the only variable in the test.Condition 0# is a blank control model.In conditions 3 #, 4 # and 5 #, the depth of karst fissure is 100 m, the inclination angle of karst fissure is 70 , and the distance from the bottom of karst fissure to the slope surface is 60 m, 90 m and 120 m, respectively.After the mining of working face, the vertical and horizontal displacement nephograms of the slope in conditions 0 #, 3#, 4#, and 5# are shown in Figures 6a, d-f and 7a, d-f, respectively.
Similarly, according to the distribution characteristics of displacement equipotential lines in the vertical and horizontal displacement nephograms, and referring to the division criteria of the slope deformed area in Section 4.1.1,the overburden above the goaf still can be divided into three parts accordance with the degree of deformation.The roof above the goaf under the four conditions all undergoes bending and subsidence, and the vertical and horizontal displacement of the trapezoidal area in the middle of the roof are about À7.94 m and À1.40 m, respectively.The displacement equipotential lines of the rock mass near the slope surface are conical distribution, and the settlement centre line of slope remains shifted in the direction of the free surface.It can be seen that the mining action is the main factor affecting the deformation of the roof above the goaf and the rock mass near the slope surface.
According to the analysis in Section 4.1.1,the existence of karst fissure directly affects the deformation of the rock mass near the slope top.As the karst fissure moves to the deep part of the slope, the vertical and horizontal deformation of the rock mass near the slope top gradually decreases.In addition to its distance from the free surface, the description of the location of karst fissure should also emphasize its relative position with the coal seam working face.Only the karst fissure located within the mining influence range is the object to be considered in the slope deformation analysis.
It can be seen more intuitively from the displacement variation curve that the vertical displacements of lines A and B across the karst fissure decrease linearly first and then decrease sharply, and the displacement remains substantially at À0.01 m (Figures 5b and 12a and b).This again shows that the karst fissure can effectively block the transmission of mining-induced deformation within the slope.The vertical deformation of lines C and D below the karst fissure first increases and then decreases from the slope surface to the inside of the slope, which also coincides with the trapezoidal or conical distribution characteristics of the displacement equipotential lines in the vertical displacement nephograms (Figures 5b and 12c and d).
In conditions 0# and 5#, the displacement amounts and change trends of the five monitoring lines are the most similar.It can be seen that the presence of the karst fissure, located in the deep part of the slope, reduce significantly the slope stability (Figure 12).

Variation law of slope stress
The vertical stress nephograms of the initial equilibrium state of the slopes and the equilibrium state of the slopes after coal seam mining in conditions 0#, 3#, 4# and 5# are shown in Figures 9a, d-f and 10a, 10d-f, respectively.
The distribution of the initial vertical stress inside the slope in conditions 4# and 5# is very similar to that in condition 3#, showing the layered distribution features of increasing compressive stress from top to bottom.The maximum compressive stress at the bottom of the model is about 7.80 Mpa.The stress equipotential lines on the left side of the bottom of the karst fissure appears obviously bulge, while the curvature of the raised stress equipotential lines is relatively flat with the increasing distance between the karst fissure and the free surface.The maximum compressive stresses at the bottom of the karst fissure in conditions 3#, 4#, and 5# are 11.70 MPa, 10.20 MPa and 9.17 MPa, respectively (Figure 9d-f).
In the mining process of working face, the mining pressurization zones are still visible on the left, middle and right side of the deformed roof.However, since the karst fissure can prevent the continuous transmission of deformation inside the slope, the stress equipotential lines within the overburden on the right side of the karst fissure are still layered.Although the position of the karst fissure is different within the slope, the position of the pressurization zone on the right side of the roof remains basically unchanged.The closer the karst fissure is to the deeper part of the slope, the larger the range of the pressurization zone on the right side of the roof is, among which the morphology of the pressurization zone on the right side of the roof is the most similar in conditions 0# and 5# (Figures 10a and 9f).By analyzing the variation law of the stress field of the slope under mining action, it can also be explained that the position of karst fissure in the slope and the relative position between karst fissure and working face are two indispensable important factors affecting the deformation of the slope.

Fracture evolution process
Taking condition 0# as the control, the influence of karst fissures at different positions in conditions 3#, 4#, and 5# on fracture evolution is analyzed with reference to Figure 11a, d-f.
After excavation of 60 m and 80 m, tension cracks and tension-shear cracks appeared for the first time at the bottom of the karst fissure and the free surface in conditions 3# and 4#, respectively (Figure 11d -2 and e -3).When excavated to 100 m, tension cracks appeared firstly at the bottom of the karst fissure in condition 5# (Figure 11f -4).After 120 m, of excavation, the tension cracks appeared at 75 m from the slope top in condition 0# (Figure 11a -5).In the mining process of working face, the deformation area of the overburden above the goaf and the fracture extension area at the bottom of the karst fissure or the top of the slope gradually approached and penetrated each other, resulting in differences in the morphology of the deformation area in the slope under different conditions (Figure 11d-f).The fracture evolution analysis also confirmed that the adverse effect of the karst fissure on the slope deformation was reduced instead when the distance of the karst fissure from the free surface was too large.respectively.The meanings of the negative values of displacement in the figures are still referred to Section 4.1.1.
Comparing the vertical and horizontal displacement nephograms of the slope in nine working conditions, it can be found that under the action of mining, except for the difference in the deformation of the rock mass on the left side of the karst fissure, the deformation features of the remaining areas in the slope are basically the same.Karst fissure, as a discontinuous plane with low strength in the rock mass, can be categorized as secondary structural plane by geological types, and also conforms to the definition of weak structural plane (Yin et al. 2021).Song (2022) pointed out that according to the relationship between the occurrence of weak structural plane and the free surface, the slope can be divided into five structural types: the directional or co-directional structure, the reverse or inverse structure, the oblique structure, the lateral structure and the flat-stack structure.The directional or co-directional structure refers to the slope structure with the same tendency of the weak structural plane and the free surface, and the strike angle of both is less than 30 ; while the reverse or inverse structure refers to the slope structure with opposite tendency of both and strike angle less than 30 .As the strike angle between the weak structural plane and the free surface gradually increases, the slope structure transitions from the oblique structure (strike angle greater than 30 and less than 60 ) to the lateral structure (strike angle greater than 60 ).When the weak structural plane is nearly horizontal to the stratum level, and the dip angle is generally not more than 10 , the slope is the flat-stack structure.Among them, sliding failure along the weak structural plane is most likely to occur on the slope with directional or co-directional structure.The above analysis also fully confirms that when the dip of the karst fissure is the same as that of the slope (70 ), the slope is the co-directional structure, and the rock mass at the slope top will slip along the long karst fissure (Figures 13a and 14a).However, with the increase of the dip of the karst fissure, the slope structure gradually changes from the directional structure to the reverse structure.The vertical displacement of rock mass at the slope top is basically maintained between À0.50 m to À2.00 m (Figure 13b-d), while the horizontal displacement gradually decreases from top to bottom (Figure 14b-d).It is inferred from the displacement variation law of rock mass that its motion state changes from sliding to toppling, and the role of karst fissure changes from a potential slip surface to the cracking boundary.
The vertical displacement variation laws of lines A to D, and the horizontal displacement variation law of line E are shown in Figure 15.The vertical displacements of the data points on the right side of the karst fissure in lines A and B are basically maintained at À0.01 m.However, except for which in condition 3# gradually decreases, the vertical displacements of the data points on the left side of the karst fissure in conditions 6 #, 7# and 8# gradually increases, and the larger the inclination of the karst fissure, the greater the vertical deformation of the data points (Figure 15a  and b).Meanwhile, the variation law of horizontal displacement of line E in conditions 6#, 7#, and 8# is also completely different from that of condition 3#, showing a gradually decreasing trend on the whole (Figure 15e).The vertical deformation of lines C and D below the karst fissure still increases first and then decreases from the slope surface to the inside of the slope, and the vertical deformation of line D is larger (Figures 5c and 15c and d).Based on the above analysis, it can be seen that the karst fissure, which has an inverse relationship with the slope surface, is also an indispensable factor to be considered in the study of slope stability.

Variation law of slope stress
The vertical stress nephograms of the initial equilibrium state of the slopes and the equilibrium state of the slopes after coal seam mining in conditions 3#, 6#, 7# and 8# are shown in Figures 16 and 17, respectively.The meanings of the positive and negative stress values are thus referred to Section 4.1.2.
As the dip of the karst fissure increases, the continuity of the initial vertical stress equipotential lines on both sides of the karst fissure gradually enhances.The dip of the karst fissure increases from 70 to 80 , 90 and 100 , the maximum compressive stress at the bottom of the karst fissure decreases abruptly from 11.70 MPa to 6.27 MPa, 2.88 MPa and 1.06 MPa, and the stress concentration zone also decreases gradually (Figure 16).The initial stress concentration zone is prone to damage under the action of mining, and is often used as the extension starting point of the mining pressurization zone on the right side of the roof.The larger the dip of the karst fissure, the smaller the range of the mining pressurization zone on the right side of the deformed roof, and the weaker the stress concentration effect (Figure 17).Therefore, it is presumed that the deformation damage of the reverse structural slope is not very serious, but still cannot be ignored.

Fracture evolution process
The fracture evolution processes of the slope in conditions 3#, 6#, 7#, and 8# during coal seam mining are shown in Figure 18.During the continuous advancement of  working face for 100 m, tension cracks appeared and developed gradually at the bottom of the karst fissure in four conditions, but tension-shear cracks appeared only at the free surface in conditions 3# and 6# (Figure 18a -3 and b -3).When a 200 m length was excavated, the deformation area of the overburden above the goaf and the fracture extension area at the bottom of the karst fissure or the top of the slope penetrated each other (Figure 18a -4 to 18d -4).According to the distribution characteristics of fissures, the deformation area of overburden above the goaf can still be divided into three parts from bottom to top: caving zone, fractured zone and bending subsidence zone.

Conclusions
Based on the modeling results and related analysis, the following conclusions are drawn.Through data collection and field investigation, it is found that most of the serious collapse and landslide disasters occurring in karst mountain areas in the southwest China are related to topography, lithology, tectonics, mining, rainfall and other factors.Pusa collapse is a typical remote chain avalanche disaster induced by underground mining.Long-term and large-scale underground mining activities accelerate the deformation and stress adjustment of the mountain.The topography and landform of steep tops and gentle bottoms provide huge potential energy for the formation of high-speed and long-distance debris flow after the instability of the mountain.The shear outlet of the leading edge of the landslide is located in the weak stratum of the lower part of the mountain.The complex karst pipeline system formed inside the mountain provides advantageous channels for rainwater infiltration.Under the action of mining, the deformation area in the slope is divided into three parts from bottom to top according to the degree of deformation: the roof of above the goaf, the rock mass near the slope surface, and the rock mass at the slope top.Fracture and collapse occurred in the middle of the roof due to excessive bending deflection.The separated fractures effectively block the upward transmission of rock strata collapse state making the vertical displacement of the rock mass near the slope surface significantly reduced.The bottom of the karst fissure is susceptible to cracking first in the process of mining due to stress concentration.And the greater the depth of the karst fissure, the more significant the phenomenon of compressive stress concentration at its bottom, and the earlier tension cracks appear.After the coal seam working face is mined, the left, middle and right sides of the deformed roof form the mining pressurization area, which coincides with the area of severe deformation in the slope.
The morphology of the karst fissure developed at the shoulder of the slope controls and determines the deformation characteristics of the rock mass at the slope top.The position of karst fissure in the slope and the relative position between karst fissure and coal seam working face are two indispensable important factors affecting the deformation of the slope.The farther the karst fissure is from the free surface, the less the adverse effect of the karst fissure on the slope deformation.Only the karst fissure located within the mining influence range is the object to be considered in the slope stability analysis.The limit karst fissure depth of the rock mass at the slope top for toppling failure is about 1/3 of the slope height.Once the depth of the karst fissure exceeds the limit depth, the karst fissure will act as a potential slip surface forcing the rock mass at the slope top to slide.With the gradual increases of the dip of the karst fissure, the relationship between the karst fissure and the free surface gradually changes from the directional or co-directional to the reverse.Meanwhile, the motion state of the rock mass at the slope top changes from sliding to toppling, and the role of karst fissure changes from a potential slip surface to the cracking boundary.Although the deformation damage of the reverse structural slope is not very serious, the influence of the karst fissure on the stability of the slope still cannot be ignored.

Disclosure statement
No potential competing interest was reported by the author(s).
Figure 1.Overview of the August 28, 2017, Pusa landslide: (a) camera image (front view); (b) camera image (rear view).Typical erosion signs at the rear edge of source area: (c) the dissolution pipelines and the dissolution fissures; (d) and (e): subsidence troughs; (f) karst cave collapse.

Figure 2 .
Figure 2. Profile of the Pusa landslide along principal sliding direction.Collated from literature Yang et al. 2020.
determined from Figure 2 that the bottom boundary elevation of the model is 1,822 m, the crest elevation is 2,147 m, and the elevations of the left and right boundaries are 1,945 m and 2,048 m, respectively.Taking the difference between the elevation of each part of the model and the elevation of the bottom boundary, the height of the top of the model slope is 325 m, and the left and right boundary heights are 123 m and 226 m, respectively.The horizontal displacement at the left and right boundaries and the vertical displacement at the bottom boundary of the slope model are constrained, and the upper boundary is the free boundary.The preliminary numerical model has 'similarity' with the actual geological conditions, and the unique geological structure of Pusa landslide is reflected in the model.The two faults set at the foot of the slope in Figure 3 are F1 and F2 in the actual slope, and four deep karst fissures at the top of the slope are developed downward.
Figure 6.Vertical displacement nephograms of the slope after the coal seam are mined for 200 m: (a) condition 0#; (b) condition 1#; (c) condition 2#; (d) condition 3#; (e) condition 4#; (f) condition 5#.The direction of the blue arrow in the figure indicates the vertical deformation direction of the roof above the goaf.

Figure 7 .
Figure 7. Horizontal displacement nephograms of the slope after the coal seam are mined for 200 m: (a) condition 0#; (b) condition 1#; (c) condition 2#; (d) condition 3#; (e) condition 4#; (f) condition 5#.The directions of the blue and red arrows in the figure indicate the horizontal deformation directions of the middle of the roof and the rock strata at the left fracture boundary, respectively.

Figure 8 .
Figure 8. Displacement variation curves of monitoring lines in slopes containing karst fissure of different depth: (a) line A; (b) line B; (c) line C; (d) line D; (e) line E.

Figure 12 .
Figure 12.Displacement variation curves of monitoring lines in slopes containing karst fissure with different position: (a) line A; (b) line B; (c) line C; (d) line D; (e) line E.

4. 3 .
Relationship between karst fissure tendency and slope deformation features4.3.1.Features of rock strata movementIn conditions 3#, 6#, 7# and 8#, the depth of karst fissure (100 m) and the distance from the bottom of the karst fissure to the slope surface (60 m) remain unchanged, and the inclination angles of karst fissure are set as 70 , 80 , 90 and 100 , respectively.When a 200 m length was excavated, the vertical and horizontal displacement nephograms of the slope in conditions 3#, 6#, 7#, and 8# are shown in Figures13 and 14 ,

Figure 13 .
Figure 13.Vertical displacement nephograms of the slope after the coal seam are mined for 200 m: (a) condition 3#; (b) condition 6#; (c) condition 7#; (d) condition 8#.The direction of the blue arrow in the figure indicates the vertical deformation direction of the roof above the goaf.

Figure 14 .
Figure 14.Horizontal displacement nephograms of the slope after the coal seam are mined for 200 m: (a) condition 3#; (b) condition 6#; (c) condition 7#; (d) condition 8#.The directions of the blue and red arrows in the figure indicate the horizontal deformation directions of the middle of the roof and the rock strata at the left fracture boundary, respectively.

Figure 15 .
Figure 15.Displacement variation curves of monitoring lines in slopes containing karst fissure with different inclination angle: (a) line A; (b) line B; (c) line C; (d) line D; (e) line E.

Table 1 .
Features of typical mining-induced landslides and collapses in karst mountain areas in the southwest China since 1980.

Table 2 .
Physico-mechanical parameters of rock mass and coal seam.

Table 3 .
Physico-mechanical parameters of rock interfaces and joints.

Table 4 .
Operation conditions of discrete element model.