# Numerical Simulation of the Effect of Stress on the Seepage of Fractured Rock Mass

Authors: Liang Shiwei; Wang Hongmei; Ge Di; Zhao Cunqing
Journal: Geofluids
 (2021)

---

## Abstract

In our country, the Yellow River Basin ecological protection and development are put forward under the background of high quality. Among the 14 large-scale coal bases in our country, 9 coal bases are located in the Yellow River Basin, and the Shenfu-Dongsheng Coalfield, which is currently the largest under development, is located here. The region is in the process of coal mining, and the movement of overlying strata will cause the stress redistribution and coal seam in overlying aquifers also due to the effect of pore water pressure along the seepage of rock fracture and damage of overlying aquifer, so in the same formation, stress and the coupled action of seepage flow will produce mutual influence. This article through the early stage of the theoretical results discussed the application of numerical simulation method for simulating 2301 face, and the effect of stress on seepage is concluded. It is proved that the numerical simulation analysis has an important reference value for the coupling problem of stress and seepage. At the same time, the protective mining of aquifers is the basic condition for the surface ecological protection of the area, and it also provides a theoretical basis for the restoration of the ecological environment in the coal mining areas of the Yellow River Basin.


---

## Body

## 1. Introduction

China’s Yellow River Basin is rich in coal resources; nine of the 14 large coal bases in our country coal base distribution are in the China Yellow River Basin’s many Dongsheng Coalfields which are located in the western Inner Mongolia, and the development of one of the biggest piece of coal is in the northern Shanxi Province, which has simple mining geological conditions; the surface has the salient feature of the fragile ecological environment. In the edge of Maowusu desert and the Loess Plateau, sparse vegetation and wind erosion are serious, and most of the areas covered by thick wind-blown sand; the main aquifer of quaternary system of Sarah WuSu group is an unconformable contact on the Jurassic bedrock; the aquifer is the survival of the whole region water system; water is very precious, mainly by the surface precipitation to supply; we have a clear understanding that the development of coal should not be at the expense of the destruction of the local ecological environment, especially in the northwest of Yulin and Ordos region; the region is an arid region, there is little rain, wind is big, and the ecological environment is very frail [1–4]. It is well known that the traditional coal mining causes damage of the underground rock mass structure, under the action of mining pressure and the movement of the overburden; the phenomenon such as fracture and overburden belongs to the discontinuous medium rock at the top of the aquifer in discontinuous medium flow, which we referred to as percolation, but due to coal mining, the influence of stress of rock mass will redistribute; it will inevitably affect the formation of the fissure zone of rock mass, which affects fluid seepage in rock mass; in the mining area, ecological restoration, protection, and restoration of aquifer are crucial. At present, scholars at home and abroad have studied more natural rock seepage, but less on the seepage characteristics under the influence of mining. In this paper, combining with the existing theory, the computer numerical simulation method commonly used in the field of geotechnical science is used to analyze this kind of problem, which has very important significance [5–8].


## 2. Brief Introduction of Previous Research Results

The relation between the permeability coefficientKh and the stress P obtained by Snow (1966) [9] is(1)Kh=K0+Aρgb24μsp−p0Kn,where K0 is the initial stress; P0 is the permeability coefficient of action; Kn is the normal stiffness of the crack; and A is the coefficient.Jones (1975) [9] proposed the empirical equation of permeability coefficient K of rock fractures as follows:(2)K=K0logphp3,where Ph is the effective stress of healing when K=0.Louis et al. (1976) [9] obtained the following according to the borehole pumping test:(3)K=K0exp−ασ,σ≈λH−p,where K0 is the surface permeability coefficient; γH is overburden weight; P is water pressure; and α is the coefficient.Kranz et al. (1979) [9] proposed the following formula to predict the change of permeability coefficient (K):(4)K∝−σc−bap.In the formula,b and a are constant.Walsh (1981) [9] obtained the empirical formula of permeability coefficient and effective stress σe as follows:(5)K=K01−2ξlnσeσeo3,where K0 is the initial permeability coefficient; σeo is the initial effective stress, and ξ is the coefficient related to the crack geometry size.Kelsall (1984) [9] proposed a method to explain the influence of stress on the permeability coefficient of rock mass:(6)Ke=Keo1+Aσeo/ξt31+Aσe/ξt3,where Ke is the permeability coefficient when the effective stress is σe and Ke acting on the σe main direction; and A, T, and ξ are the material constants of rock mass fractures.Bai and Elsworth (1989) [9] assumed that cracks were soft, and the following equation was determined:(7)ΔK=ρg12sμb+sΔε3,where Δε is the strain perpendicular to the fracture group.Zhang et al. [10] obtained the relationship between pressure, permeability coefficient, and seepage flow as follows:(8)K=K01+p−γHcos2θ+λsin2θbKn4,Q=Q01+p−γHcos2θ+λsin2θbKn4.


## 3. Stress Seepage Coupling Numerical Simulation Analysis of Fractured Rock Mass

### 3.1. Establishment of Mechanical Model

The 2301 working face of a mine is located in the desert area of Yulin area, the ecological environment is fragile, the working face is located in the 3 coal seam, the average thickness is 3.5 meters, and the overburden contains aquifer, which has an important impact on the mining of coal mine; improper mining will cause the accident of water inrush and sand. According to the rock mechanics experiment, the rock mechanics parameters are obtained (Table1).Table 1
The face of 2301 of rock physical and mechanical properties.
Serial numberRock nameLayer thickness (m)Buried depth (m)Bulk density (kg/m3)Compressive strength (MPa)Tensile strength (MPa)Elastic modulus (GPa)1Weathered mudstone4343230010.61.21.02Sandstone104.2147.2255031.31.27.63Sand and mudstone1.7148.9273066.71.86.243 coal seams3.5152.41420270.45Sandy mudstone40.2192.6261035.32.81.2According to the roadway layout characteristics of 2301 working face and the near-horizontal coal seam (as shown in Figure1), in order to improve the calculation efficiency, the design model is shown in Figure 2. This simulation is mainly to analyze the distribution law and development height of caving zone and fracture zone with overlying aquifers. The entire overlying rock and a small part of the floor should be considered when designing the model, which is to analyze the coordination of coal seam, roof, and floor; this is key to analyze the problem. This thesis has certain pore water pressure of the aquifer as module, located 10 m beneath the floor of the overlying strata in the process of working face advancing; fluid-structure interaction phenomenon has mutual influence, the change of pore water pressure in the rock mass causes stress variation in the mining, and water may exist in the crack, and the crack volume changes; it will have an effect on pore water pressure; this model uses the fluid calculation and mechanical calculation; at the same time, it also is often said that the fluid-structure coupling calculation adopts the model of seepage flow in the calculation and isotropic seepage model, calculating result in order to achieve. In the model, the range of coal seam roof above is until 148.9 m as a model of the upper surface and a layer of rock as a model of the base plate; each side of the model scales out 50 m; boundary treatment is the key problem of the model; the specific process is as follows: on the left and right boundary of the model, take u=0, v=0 (u for the X axis displacement, v for the Y direction displacement, namely, single constraint boundary). At the bottom boundary of the model, u=v=0 is the fully constrained boundary. The upper boundary of the model is the ground surface [11–15], which is regarded as free boundary without constraint [16–18].Figure 1
Roadway layout of working face.Figure 2
The establishment of the mechanical model. 1 sandy mudstone; 2 coal seams; 3 sand and mudstone; 4 sandstone; 5 weathered mudstone.


### 3.2. Coupling Analysis of Stress and Seepage

The movement of the overlying strata in coal mining leads to the redistribution of stress, which is bound to affect the flow of fluid in the fracture, and the key parameter is the permeability coefficient, which reflects the flow law of fluid and also reflects the development law of fracture in the overlying strata.


#### 3.2.1. Horizontal Movement Law of Overburden

Due to the influence of working face mining, each unit of the mathematical model has to undergo continuous deformation. In order to clearly display the horizontal displacement change process of each node unit, the displacement calculation results of each node were output along with the working face advancing, and the displacement cloud maps of 2301 working face with different lengths and different stoping stages were drawn. The horizontal displacement cloud maps of working face advancing with 100 m, 150 m, 200 m, and 250 m are shown in Figures3–5 respectively.Figure 3
Horizontal displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m horizontal displacement cloud map(b)
150 m horizontal displacement cloud map(c)
200 m horizontal displacement cloud map(d)
250 m horizontal displacement cloud mapFigure 4
Horizontal displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 220 m.
(a)
100 m horizontal displacement cloud map(b)
150 m horizontal displacement cloud map(c)
200 m horizontal displacement cloud map(d)
250 m horizontal displacement cloud mapFigure 5
Horizontal displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 240 m.
(a)
100 m horizontal displacement cloud map(b)
150 m horizontal displacement cloud map(c)
200 m horizontal displacement cloud map(d)
250 m horizontal displacement cloud mapIn the working face, the phase of initial overall overburden deformation is not large, including horizontal deformation value; this is due to the fact that the immediate roof is not fully caving; the supporting role, on both sides of the mined-out area, has the opposite change; there is displacement of the change in the opposite direction; its value is nearer the mined-out area; the greater the distance of the far mined-out area, the more and more small, along with the mining, mine pressure release, the overburden, and the horizontal displacement of zero, sometimes at the same time, due to the restriction of the boundary conditions, the model, and the lower part of the border around the border; the strata horizontal displacement value is separated into four areas: there are negative moving area, positive moving area far away from the goaf, and positive moving area and negative moving area near the goaf.


#### 3.2.2. Vertical Movement Law of Overburden

In the different stages of 2301 working face, the vertical subsidence curve shape of overlying rocks does not change much in the process of working face advancement, and the changes of each rock layer are basically the same; whether it is hard rock or soft rock, all are coordinated deformation.In the overburden on the vertical displacement of the funnel characteristics, in early mining, the vertical displacement is not big, but the amount of propulsion, strata of vertical displacement, and the maximum value appeared in the middle of the mined-out area. This is because in the middle of the mined-out area, without support, the rock strata have a greater range of movement, and from the picture that looks like a funnel, when mining fully, strata are not a big change, its change of displacement distribution is symmetrical, goaf in the working face advances continuously in the process of filling, and vertical displacement of overburden rock will then move forward, until it is focused on the coal wall.From the vertical displacement cloud maps of different mining stages, in the beginning stage, the displacement cloud maps in the vertical direction are not very large either in size or in advance influence range. With the increase of distance, these changes also increase, and the vertical displacement cloud map gradually becomes unclosed (Figures6–8).Figure 6
Vertical displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m vertical displacement cloud map(b)
150 m vertical displacement cloud map(c)
200 m vertical displacement cloud map(d)
250 m vertical displacement cloud mapFigure 7
Vertical displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 220 m.
(a)
100 m vertical displacement cloud map(b)
150 m vertical displacement cloud map(c)
200 m vertical displacement cloud map(d)
250 m vertical displacement cloud mapFigure 8
Vertical displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 240 m.
(a)
100 m vertical displacement cloud map(b)
150 m vertical displacement cloud map(c)
200 m vertical displacement cloud map(d)
250 m vertical displacement cloud map


#### 3.2.3. Horizontal Stress Distribution in Rock Mass

Numerical simulation results show that the early mining overburden horizontal stress is not big. However, the zoning of the overlying strata in the whole model is obvious. This is because the overlying strata are composed of rock layers of different properties during the deposition process. At the same time, the vertical distance between each layer of rock strata and the working surface is also different, so in the face of the disturbance of underground mining, it will show different characteristics of overlying rock changes. The distribution of strata stress is stress arch, usually in place after the coal for starting cut arch foot, in front of the coal wall stress peak area of stress arch structure of the former arch foot; with the working face advancing, the arch is not fixed but constantly moves forward.After the end of mining, the horizontal stress distribution in the mined-out area is as follows: the horizontal stress in the mined-out area is larger 20~40 m upward, indicating that the influence of stress mining is greater and strong. The change of horizontal stress in the goaf upward 40-80 m is moderate, and the fluctuation is not large. The horizontal stress from above 80 m to the surface is small. It is due to the mining of the underground working face, which is far away from the surface of the ground, and the direct roof and basic roof rock layers below have undertaken most of the horizontal deformation. As you can see from the following diagram, there is a significant partitioning phenomenon. It also confirms that the extent of damage to the overlying rock strata caused by mining activities gradually weakens with the increase of distance (Figures9–11).Figure 9
Horizontal stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m horizontal stress nephogram(b)
150 m horizontal stress nephogram(c)
200 m horizontal stress nephogram(d)
250 m horizontal stress nephogramFigure 10
Horizontal stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 220 m.
(a)
100 m horizontal stress nephogram(b)
150 m horizontal stress nephogram(c)
200 m horizontal stress nephogram(d)
250 m horizontal stress nephogramFigure 11
Horizontal stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 240 m.
(a)
100 m horizontal stress nephogram(b)
150 m horizon al stress nephogram(c)
200 m horizontal stress nephogram(d)
250 m horizontal stress nephogram


#### 3.2.4. Distribution Law of Vertical Stress inside Rock Mass

In the process of advancing the working face, the overburden undergoes the process of “deformation-separation-instability.” At the same time, the internal stress of the rock mass also changes accordingly. Because the basic roof is thick and hard and the ability to resist deformation is strong, the duration of each stage is longer. When the working face is advanced to 150 m, the advance influence range of the working face abutment pressure is about 50 m. With the advance of working face, the influence range of leading pressure is also changing, and it has experienced the process of “increase-decrease-increase.” After mining, stress concentration appears especially near the top of the coal pillar. The peak depth of the coal wall is 10 m, and the influence range of abutment pressure is about 60 m in front of the coal wall (Figures12–14).Figure 12
Vertical stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m vertical stress nephogram(b)
150 m vertical stress nephogram(c)
200 m vertical stress nephogram(d)
250 m vertical stress nephogramFigure 13
Vertical stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 220 m.
(a)
100 m vertical stress nephogram(b)
150 m vertical stress nephogram(c)
200 m vertical stress nephogram(d)
250 m vertical stress nephogramFigure 14
Vertical stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 240 m.
(a)
100 m vertical stress nephogram(b)
150 m vertical stress nephogram(c)
200 m vertical stress nephogram(d)
250 m vertical stress nephogram


#### 3.2.5. Overburden Failure Status and Plastic Zone Distribution

In the elastoplastic analysis, when the FLAC program is dealing with the problem, it will judge the stress state of the point according to the Mohr-Coulomb criterion at the same time when calculating the stress and displacement of each element. In the numerical simulation, the following conventions are made: 0 denoted that the rock mass element is still in elastic state;1: indicates that the rock mass element is yielding and begins to transition to the plastic state; 2 represents that the rock mass element enters the plastic state; 3 means that the rock mass element has exceeded the uniaxial tensile resistance; 4 represents the yield of rock mass unit and exceeds uniaxial tensile resistance; and 5 indicates that the rock mass element has exceeded the uniaxial tensile resistance. The software uses the existing failure parameters of each unit to integrate, so as to express the failure situation of the whole overburden. The identification of plastic zone is the main basis of the failure state of overburden, as well as the method to judge the distribution law and development height of the water-conducting fracture zone.It can be seen from the plastic separation diagram that the main failure places are at the top and left and right sides of the model. Under the influence of mining damage, the overburden failure is symmetrically distributed. By comprehensive comparison of various calculation models, different working face lengths lead to different distribution of plastic zone and development height of the water-conducting fracture zone. According to the analysis in Figure15, when the working face length is 200 m, the maximum height of plastic zone distribution in overlying rock is 48 m. When the working face length of 220 m strata distribution of plastic zone is in the maximum height of 54 m and when the working face length of 240 m strata distribution of plastic zone in the maximum height has communicated the overlying loose bed, the overlying aquifer water inrush occurs; in production practice, it cannot meet the working face length. Under the action of mining seepage, the length of the working face cannot be expanded indefinitely. Once it exceeds a certain value, the overburden fissures will communicate with the aquifer and even the surface, causing water inrush accidents. So the working face length of 220 m is the most reasonable [19–22].Figure 15
The plastic zone distribution when the working face length is 200 m, 220 m, 240 m, and pushed to 250 m.
(a)
Working face length is 200 m and advances to 250 m(b)
Working face length is 220 m to 250 m(c)
Working face length is 240 m advancing to 250 m


#### 3.2.6. Coupling Analysis of Stress and Seepage of Overlying Rock when Working Face Length is 220 m

Under the influence of mining, the stress of the overlying strata changes, which makes the permeability of the overlying strata change. Therefore, the coupling analysis of stress and seepage is particularly important in engineering applications. The seepage problem of overburden rock caused by mining is mainly to study the law of fluid movement in pores and fissures. Among them, the parameter that best reflects the law of change is the permeability coefficient. The change in permeability coefficient reflects the expansion and closure of pores and fissures. This also explains the overburden from the side. In the law of rock movement, we can analyze the law of seepage through the change of stress, that is, the coupling analysis of stress and seepage.Here, we take the length of the working face 220 m as the basis. In the mining geological conditions, we know that there is an aquifer above the bedrock, which will inevitably generate pore water pressure on the surrounding rock formations. The existence of pore water pressure will affect the surrounding rocks. The nature of the pore water will change, and the change of pore water pressure will affect the stress and strain of the rock formation, which will redistribute the stress of the rock formation. In the mining of the working face, the upper and lower ends of the working face and the cut-off are most affected. At this time, the permeability of the rock formation will increase accordingly. As the work surface advances, the stress of the overlying strata of the stope will be redistributed, and the overlying strata will be closed from the original tensile cracks and the permeability at this time The rate gradually decreases until the fissure is closed, the pore water pressure is zero, the working face continues to advance, and the pore water pressure may become negative, because the water-conducting fissure zone may spread to the overlying aquifer, affecting the water level of the aquifer.According to the simulation, the distribution of pore water pressure reflects the movement of the overlying strata in the process of advancing the working face. The figure shows that the stress distribution at both ends of the working face is concentrated, and a falling funnel appears above the goaf (Figure16).Figure 16
Cloud diagram of overburden pore water pressure at 250 m advance of working face advances to 100 m, 150 m, 200 m, and 250 m.
(a)
100 m advance of working face(b)
150 m advance of working face(c)
200 m advance of working face(d)
250 m advance of working faceAccording to the numerical simulation results, the change rule of overburden stress under the action of overburden pore water pressure can be obtained, which is expressed by the trend in Figure17.Figure 17
Comparison of permeability coefficient variation trend of different advancing distances at different working faces.
(a)(b)(c)Through the above analysis, the movement fracture of the overlying strata in coal mining is the result of the joint action of stress and seepage. The change of stress will change the flow law of fluid accordingly.


## 3.1. Establishment of Mechanical Model

The 2301 working face of a mine is located in the desert area of Yulin area, the ecological environment is fragile, the working face is located in the 3 coal seam, the average thickness is 3.5 meters, and the overburden contains aquifer, which has an important impact on the mining of coal mine; improper mining will cause the accident of water inrush and sand. According to the rock mechanics experiment, the rock mechanics parameters are obtained (Table1).Table 1
The face of 2301 of rock physical and mechanical properties.
Serial numberRock nameLayer thickness (m)Buried depth (m)Bulk density (kg/m3)Compressive strength (MPa)Tensile strength (MPa)Elastic modulus (GPa)1Weathered mudstone4343230010.61.21.02Sandstone104.2147.2255031.31.27.63Sand and mudstone1.7148.9273066.71.86.243 coal seams3.5152.41420270.45Sandy mudstone40.2192.6261035.32.81.2According to the roadway layout characteristics of 2301 working face and the near-horizontal coal seam (as shown in Figure1), in order to improve the calculation efficiency, the design model is shown in Figure 2. This simulation is mainly to analyze the distribution law and development height of caving zone and fracture zone with overlying aquifers. The entire overlying rock and a small part of the floor should be considered when designing the model, which is to analyze the coordination of coal seam, roof, and floor; this is key to analyze the problem. This thesis has certain pore water pressure of the aquifer as module, located 10 m beneath the floor of the overlying strata in the process of working face advancing; fluid-structure interaction phenomenon has mutual influence, the change of pore water pressure in the rock mass causes stress variation in the mining, and water may exist in the crack, and the crack volume changes; it will have an effect on pore water pressure; this model uses the fluid calculation and mechanical calculation; at the same time, it also is often said that the fluid-structure coupling calculation adopts the model of seepage flow in the calculation and isotropic seepage model, calculating result in order to achieve. In the model, the range of coal seam roof above is until 148.9 m as a model of the upper surface and a layer of rock as a model of the base plate; each side of the model scales out 50 m; boundary treatment is the key problem of the model; the specific process is as follows: on the left and right boundary of the model, take u=0, v=0 (u for the X axis displacement, v for the Y direction displacement, namely, single constraint boundary). At the bottom boundary of the model, u=v=0 is the fully constrained boundary. The upper boundary of the model is the ground surface [11–15], which is regarded as free boundary without constraint [16–18].Figure 1
Roadway layout of working face.Figure 2
The establishment of the mechanical model. 1 sandy mudstone; 2 coal seams; 3 sand and mudstone; 4 sandstone; 5 weathered mudstone.


## 3.2. Coupling Analysis of Stress and Seepage

The movement of the overlying strata in coal mining leads to the redistribution of stress, which is bound to affect the flow of fluid in the fracture, and the key parameter is the permeability coefficient, which reflects the flow law of fluid and also reflects the development law of fracture in the overlying strata.


### 3.2.1. Horizontal Movement Law of Overburden

Due to the influence of working face mining, each unit of the mathematical model has to undergo continuous deformation. In order to clearly display the horizontal displacement change process of each node unit, the displacement calculation results of each node were output along with the working face advancing, and the displacement cloud maps of 2301 working face with different lengths and different stoping stages were drawn. The horizontal displacement cloud maps of working face advancing with 100 m, 150 m, 200 m, and 250 m are shown in Figures3–5 respectively.Figure 3
Horizontal displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m horizontal displacement cloud map(b)
150 m horizontal displacement cloud map(c)
200 m horizontal displacement cloud map(d)
250 m horizontal displacement cloud mapFigure 4
Horizontal displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 220 m.
(a)
100 m horizontal displacement cloud map(b)
150 m horizontal displacement cloud map(c)
200 m horizontal displacement cloud map(d)
250 m horizontal displacement cloud mapFigure 5
Horizontal displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 240 m.
(a)
100 m horizontal displacement cloud map(b)
150 m horizontal displacement cloud map(c)
200 m horizontal displacement cloud map(d)
250 m horizontal displacement cloud mapIn the working face, the phase of initial overall overburden deformation is not large, including horizontal deformation value; this is due to the fact that the immediate roof is not fully caving; the supporting role, on both sides of the mined-out area, has the opposite change; there is displacement of the change in the opposite direction; its value is nearer the mined-out area; the greater the distance of the far mined-out area, the more and more small, along with the mining, mine pressure release, the overburden, and the horizontal displacement of zero, sometimes at the same time, due to the restriction of the boundary conditions, the model, and the lower part of the border around the border; the strata horizontal displacement value is separated into four areas: there are negative moving area, positive moving area far away from the goaf, and positive moving area and negative moving area near the goaf.


### 3.2.2. Vertical Movement Law of Overburden

In the different stages of 2301 working face, the vertical subsidence curve shape of overlying rocks does not change much in the process of working face advancement, and the changes of each rock layer are basically the same; whether it is hard rock or soft rock, all are coordinated deformation.In the overburden on the vertical displacement of the funnel characteristics, in early mining, the vertical displacement is not big, but the amount of propulsion, strata of vertical displacement, and the maximum value appeared in the middle of the mined-out area. This is because in the middle of the mined-out area, without support, the rock strata have a greater range of movement, and from the picture that looks like a funnel, when mining fully, strata are not a big change, its change of displacement distribution is symmetrical, goaf in the working face advances continuously in the process of filling, and vertical displacement of overburden rock will then move forward, until it is focused on the coal wall.From the vertical displacement cloud maps of different mining stages, in the beginning stage, the displacement cloud maps in the vertical direction are not very large either in size or in advance influence range. With the increase of distance, these changes also increase, and the vertical displacement cloud map gradually becomes unclosed (Figures6–8).Figure 6
Vertical displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m vertical displacement cloud map(b)
150 m vertical displacement cloud map(c)
200 m vertical displacement cloud map(d)
250 m vertical displacement cloud mapFigure 7
Vertical displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 220 m.
(a)
100 m vertical displacement cloud map(b)
150 m vertical displacement cloud map(c)
200 m vertical displacement cloud map(d)
250 m vertical displacement cloud mapFigure 8
Vertical displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 240 m.
(a)
100 m vertical displacement cloud map(b)
150 m vertical displacement cloud map(c)
200 m vertical displacement cloud map(d)
250 m vertical displacement cloud map


### 3.2.3. Horizontal Stress Distribution in Rock Mass

Numerical simulation results show that the early mining overburden horizontal stress is not big. However, the zoning of the overlying strata in the whole model is obvious. This is because the overlying strata are composed of rock layers of different properties during the deposition process. At the same time, the vertical distance between each layer of rock strata and the working surface is also different, so in the face of the disturbance of underground mining, it will show different characteristics of overlying rock changes. The distribution of strata stress is stress arch, usually in place after the coal for starting cut arch foot, in front of the coal wall stress peak area of stress arch structure of the former arch foot; with the working face advancing, the arch is not fixed but constantly moves forward.After the end of mining, the horizontal stress distribution in the mined-out area is as follows: the horizontal stress in the mined-out area is larger 20~40 m upward, indicating that the influence of stress mining is greater and strong. The change of horizontal stress in the goaf upward 40-80 m is moderate, and the fluctuation is not large. The horizontal stress from above 80 m to the surface is small. It is due to the mining of the underground working face, which is far away from the surface of the ground, and the direct roof and basic roof rock layers below have undertaken most of the horizontal deformation. As you can see from the following diagram, there is a significant partitioning phenomenon. It also confirms that the extent of damage to the overlying rock strata caused by mining activities gradually weakens with the increase of distance (Figures9–11).Figure 9
Horizontal stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m horizontal stress nephogram(b)
150 m horizontal stress nephogram(c)
200 m horizontal stress nephogram(d)
250 m horizontal stress nephogramFigure 10
Horizontal stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 220 m.
(a)
100 m horizontal stress nephogram(b)
150 m horizontal stress nephogram(c)
200 m horizontal stress nephogram(d)
250 m horizontal stress nephogramFigure 11
Horizontal stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 240 m.
(a)
100 m horizontal stress nephogram(b)
150 m horizon al stress nephogram(c)
200 m horizontal stress nephogram(d)
250 m horizontal stress nephogram


### 3.2.4. Distribution Law of Vertical Stress inside Rock Mass

In the process of advancing the working face, the overburden undergoes the process of “deformation-separation-instability.” At the same time, the internal stress of the rock mass also changes accordingly. Because the basic roof is thick and hard and the ability to resist deformation is strong, the duration of each stage is longer. When the working face is advanced to 150 m, the advance influence range of the working face abutment pressure is about 50 m. With the advance of working face, the influence range of leading pressure is also changing, and it has experienced the process of “increase-decrease-increase.” After mining, stress concentration appears especially near the top of the coal pillar. The peak depth of the coal wall is 10 m, and the influence range of abutment pressure is about 60 m in front of the coal wall (Figures12–14).Figure 12
Vertical stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m vertical stress nephogram(b)
150 m vertical stress nephogram(c)
200 m vertical stress nephogram(d)
250 m vertical stress nephogramFigure 13
Vertical stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 220 m.
(a)
100 m vertical stress nephogram(b)
150 m vertical stress nephogram(c)
200 m vertical stress nephogram(d)
250 m vertical stress nephogramFigure 14
Vertical stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 240 m.
(a)
100 m vertical stress nephogram(b)
150 m vertical stress nephogram(c)
200 m vertical stress nephogram(d)
250 m vertical stress nephogram


### 3.2.5. Overburden Failure Status and Plastic Zone Distribution

In the elastoplastic analysis, when the FLAC program is dealing with the problem, it will judge the stress state of the point according to the Mohr-Coulomb criterion at the same time when calculating the stress and displacement of each element. In the numerical simulation, the following conventions are made: 0 denoted that the rock mass element is still in elastic state;1: indicates that the rock mass element is yielding and begins to transition to the plastic state; 2 represents that the rock mass element enters the plastic state; 3 means that the rock mass element has exceeded the uniaxial tensile resistance; 4 represents the yield of rock mass unit and exceeds uniaxial tensile resistance; and 5 indicates that the rock mass element has exceeded the uniaxial tensile resistance. The software uses the existing failure parameters of each unit to integrate, so as to express the failure situation of the whole overburden. The identification of plastic zone is the main basis of the failure state of overburden, as well as the method to judge the distribution law and development height of the water-conducting fracture zone.It can be seen from the plastic separation diagram that the main failure places are at the top and left and right sides of the model. Under the influence of mining damage, the overburden failure is symmetrically distributed. By comprehensive comparison of various calculation models, different working face lengths lead to different distribution of plastic zone and development height of the water-conducting fracture zone. According to the analysis in Figure15, when the working face length is 200 m, the maximum height of plastic zone distribution in overlying rock is 48 m. When the working face length of 220 m strata distribution of plastic zone is in the maximum height of 54 m and when the working face length of 240 m strata distribution of plastic zone in the maximum height has communicated the overlying loose bed, the overlying aquifer water inrush occurs; in production practice, it cannot meet the working face length. Under the action of mining seepage, the length of the working face cannot be expanded indefinitely. Once it exceeds a certain value, the overburden fissures will communicate with the aquifer and even the surface, causing water inrush accidents. So the working face length of 220 m is the most reasonable [19–22].Figure 15
The plastic zone distribution when the working face length is 200 m, 220 m, 240 m, and pushed to 250 m.
(a)
Working face length is 200 m and advances to 250 m(b)
Working face length is 220 m to 250 m(c)
Working face length is 240 m advancing to 250 m


### 3.2.6. Coupling Analysis of Stress and Seepage of Overlying Rock when Working Face Length is 220 m

Under the influence of mining, the stress of the overlying strata changes, which makes the permeability of the overlying strata change. Therefore, the coupling analysis of stress and seepage is particularly important in engineering applications. The seepage problem of overburden rock caused by mining is mainly to study the law of fluid movement in pores and fissures. Among them, the parameter that best reflects the law of change is the permeability coefficient. The change in permeability coefficient reflects the expansion and closure of pores and fissures. This also explains the overburden from the side. In the law of rock movement, we can analyze the law of seepage through the change of stress, that is, the coupling analysis of stress and seepage.Here, we take the length of the working face 220 m as the basis. In the mining geological conditions, we know that there is an aquifer above the bedrock, which will inevitably generate pore water pressure on the surrounding rock formations. The existence of pore water pressure will affect the surrounding rocks. The nature of the pore water will change, and the change of pore water pressure will affect the stress and strain of the rock formation, which will redistribute the stress of the rock formation. In the mining of the working face, the upper and lower ends of the working face and the cut-off are most affected. At this time, the permeability of the rock formation will increase accordingly. As the work surface advances, the stress of the overlying strata of the stope will be redistributed, and the overlying strata will be closed from the original tensile cracks and the permeability at this time The rate gradually decreases until the fissure is closed, the pore water pressure is zero, the working face continues to advance, and the pore water pressure may become negative, because the water-conducting fissure zone may spread to the overlying aquifer, affecting the water level of the aquifer.According to the simulation, the distribution of pore water pressure reflects the movement of the overlying strata in the process of advancing the working face. The figure shows that the stress distribution at both ends of the working face is concentrated, and a falling funnel appears above the goaf (Figure16).Figure 16
Cloud diagram of overburden pore water pressure at 250 m advance of working face advances to 100 m, 150 m, 200 m, and 250 m.
(a)
100 m advance of working face(b)
150 m advance of working face(c)
200 m advance of working face(d)
250 m advance of working faceAccording to the numerical simulation results, the change rule of overburden stress under the action of overburden pore water pressure can be obtained, which is expressed by the trend in Figure17.Figure 17
Comparison of permeability coefficient variation trend of different advancing distances at different working faces.
(a)(b)(c)Through the above analysis, the movement fracture of the overlying strata in coal mining is the result of the joint action of stress and seepage. The change of stress will change the flow law of fluid accordingly.


## 3.2.1. Horizontal Movement Law of Overburden

Due to the influence of working face mining, each unit of the mathematical model has to undergo continuous deformation. In order to clearly display the horizontal displacement change process of each node unit, the displacement calculation results of each node were output along with the working face advancing, and the displacement cloud maps of 2301 working face with different lengths and different stoping stages were drawn. The horizontal displacement cloud maps of working face advancing with 100 m, 150 m, 200 m, and 250 m are shown in Figures3–5 respectively.Figure 3
Horizontal displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m horizontal displacement cloud map(b)
150 m horizontal displacement cloud map(c)
200 m horizontal displacement cloud map(d)
250 m horizontal displacement cloud mapFigure 4
Horizontal displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 220 m.
(a)
100 m horizontal displacement cloud map(b)
150 m horizontal displacement cloud map(c)
200 m horizontal displacement cloud map(d)
250 m horizontal displacement cloud mapFigure 5
Horizontal displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 240 m.
(a)
100 m horizontal displacement cloud map(b)
150 m horizontal displacement cloud map(c)
200 m horizontal displacement cloud map(d)
250 m horizontal displacement cloud mapIn the working face, the phase of initial overall overburden deformation is not large, including horizontal deformation value; this is due to the fact that the immediate roof is not fully caving; the supporting role, on both sides of the mined-out area, has the opposite change; there is displacement of the change in the opposite direction; its value is nearer the mined-out area; the greater the distance of the far mined-out area, the more and more small, along with the mining, mine pressure release, the overburden, and the horizontal displacement of zero, sometimes at the same time, due to the restriction of the boundary conditions, the model, and the lower part of the border around the border; the strata horizontal displacement value is separated into four areas: there are negative moving area, positive moving area far away from the goaf, and positive moving area and negative moving area near the goaf.


## 3.2.2. Vertical Movement Law of Overburden

In the different stages of 2301 working face, the vertical subsidence curve shape of overlying rocks does not change much in the process of working face advancement, and the changes of each rock layer are basically the same; whether it is hard rock or soft rock, all are coordinated deformation.In the overburden on the vertical displacement of the funnel characteristics, in early mining, the vertical displacement is not big, but the amount of propulsion, strata of vertical displacement, and the maximum value appeared in the middle of the mined-out area. This is because in the middle of the mined-out area, without support, the rock strata have a greater range of movement, and from the picture that looks like a funnel, when mining fully, strata are not a big change, its change of displacement distribution is symmetrical, goaf in the working face advances continuously in the process of filling, and vertical displacement of overburden rock will then move forward, until it is focused on the coal wall.From the vertical displacement cloud maps of different mining stages, in the beginning stage, the displacement cloud maps in the vertical direction are not very large either in size or in advance influence range. With the increase of distance, these changes also increase, and the vertical displacement cloud map gradually becomes unclosed (Figures6–8).Figure 6
Vertical displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m vertical displacement cloud map(b)
150 m vertical displacement cloud map(c)
200 m vertical displacement cloud map(d)
250 m vertical displacement cloud mapFigure 7
Vertical displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 220 m.
(a)
100 m vertical displacement cloud map(b)
150 m vertical displacement cloud map(c)
200 m vertical displacement cloud map(d)
250 m vertical displacement cloud mapFigure 8
Vertical displacement cloud map of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 240 m.
(a)
100 m vertical displacement cloud map(b)
150 m vertical displacement cloud map(c)
200 m vertical displacement cloud map(d)
250 m vertical displacement cloud map


## 3.2.3. Horizontal Stress Distribution in Rock Mass

Numerical simulation results show that the early mining overburden horizontal stress is not big. However, the zoning of the overlying strata in the whole model is obvious. This is because the overlying strata are composed of rock layers of different properties during the deposition process. At the same time, the vertical distance between each layer of rock strata and the working surface is also different, so in the face of the disturbance of underground mining, it will show different characteristics of overlying rock changes. The distribution of strata stress is stress arch, usually in place after the coal for starting cut arch foot, in front of the coal wall stress peak area of stress arch structure of the former arch foot; with the working face advancing, the arch is not fixed but constantly moves forward.After the end of mining, the horizontal stress distribution in the mined-out area is as follows: the horizontal stress in the mined-out area is larger 20~40 m upward, indicating that the influence of stress mining is greater and strong. The change of horizontal stress in the goaf upward 40-80 m is moderate, and the fluctuation is not large. The horizontal stress from above 80 m to the surface is small. It is due to the mining of the underground working face, which is far away from the surface of the ground, and the direct roof and basic roof rock layers below have undertaken most of the horizontal deformation. As you can see from the following diagram, there is a significant partitioning phenomenon. It also confirms that the extent of damage to the overlying rock strata caused by mining activities gradually weakens with the increase of distance (Figures9–11).Figure 9
Horizontal stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m horizontal stress nephogram(b)
150 m horizontal stress nephogram(c)
200 m horizontal stress nephogram(d)
250 m horizontal stress nephogramFigure 10
Horizontal stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 220 m.
(a)
100 m horizontal stress nephogram(b)
150 m horizontal stress nephogram(c)
200 m horizontal stress nephogram(d)
250 m horizontal stress nephogramFigure 11
Horizontal stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the working face length is 240 m.
(a)
100 m horizontal stress nephogram(b)
150 m horizon al stress nephogram(c)
200 m horizontal stress nephogram(d)
250 m horizontal stress nephogram


## 3.2.4. Distribution Law of Vertical Stress inside Rock Mass

In the process of advancing the working face, the overburden undergoes the process of “deformation-separation-instability.” At the same time, the internal stress of the rock mass also changes accordingly. Because the basic roof is thick and hard and the ability to resist deformation is strong, the duration of each stage is longer. When the working face is advanced to 150 m, the advance influence range of the working face abutment pressure is about 50 m. With the advance of working face, the influence range of leading pressure is also changing, and it has experienced the process of “increase-decrease-increase.” After mining, stress concentration appears especially near the top of the coal pillar. The peak depth of the coal wall is 10 m, and the influence range of abutment pressure is about 60 m in front of the coal wall (Figures12–14).Figure 12
Vertical stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when the length of working face is 200 m.
(a)
100 m vertical stress nephogram(b)
150 m vertical stress nephogram(c)
200 m vertical stress nephogram(d)
250 m vertical stress nephogramFigure 13
Vertical stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 220 m.
(a)
100 m vertical stress nephogram(b)
150 m vertical stress nephogram(c)
200 m vertical stress nephogram(d)
250 m vertical stress nephogramFigure 14
Vertical stress nephogram of different stoping stages (100 m, 150 m, 200 m, and 250 m) when working face length is 240 m.
(a)
100 m vertical stress nephogram(b)
150 m vertical stress nephogram(c)
200 m vertical stress nephogram(d)
250 m vertical stress nephogram


## 3.2.5. Overburden Failure Status and Plastic Zone Distribution

In the elastoplastic analysis, when the FLAC program is dealing with the problem, it will judge the stress state of the point according to the Mohr-Coulomb criterion at the same time when calculating the stress and displacement of each element. In the numerical simulation, the following conventions are made: 0 denoted that the rock mass element is still in elastic state;1: indicates that the rock mass element is yielding and begins to transition to the plastic state; 2 represents that the rock mass element enters the plastic state; 3 means that the rock mass element has exceeded the uniaxial tensile resistance; 4 represents the yield of rock mass unit and exceeds uniaxial tensile resistance; and 5 indicates that the rock mass element has exceeded the uniaxial tensile resistance. The software uses the existing failure parameters of each unit to integrate, so as to express the failure situation of the whole overburden. The identification of plastic zone is the main basis of the failure state of overburden, as well as the method to judge the distribution law and development height of the water-conducting fracture zone.It can be seen from the plastic separation diagram that the main failure places are at the top and left and right sides of the model. Under the influence of mining damage, the overburden failure is symmetrically distributed. By comprehensive comparison of various calculation models, different working face lengths lead to different distribution of plastic zone and development height of the water-conducting fracture zone. According to the analysis in Figure15, when the working face length is 200 m, the maximum height of plastic zone distribution in overlying rock is 48 m. When the working face length of 220 m strata distribution of plastic zone is in the maximum height of 54 m and when the working face length of 240 m strata distribution of plastic zone in the maximum height has communicated the overlying loose bed, the overlying aquifer water inrush occurs; in production practice, it cannot meet the working face length. Under the action of mining seepage, the length of the working face cannot be expanded indefinitely. Once it exceeds a certain value, the overburden fissures will communicate with the aquifer and even the surface, causing water inrush accidents. So the working face length of 220 m is the most reasonable [19–22].Figure 15
The plastic zone distribution when the working face length is 200 m, 220 m, 240 m, and pushed to 250 m.
(a)
Working face length is 200 m and advances to 250 m(b)
Working face length is 220 m to 250 m(c)
Working face length is 240 m advancing to 250 m


## 3.2.6. Coupling Analysis of Stress and Seepage of Overlying Rock when Working Face Length is 220 m

Under the influence of mining, the stress of the overlying strata changes, which makes the permeability of the overlying strata change. Therefore, the coupling analysis of stress and seepage is particularly important in engineering applications. The seepage problem of overburden rock caused by mining is mainly to study the law of fluid movement in pores and fissures. Among them, the parameter that best reflects the law of change is the permeability coefficient. The change in permeability coefficient reflects the expansion and closure of pores and fissures. This also explains the overburden from the side. In the law of rock movement, we can analyze the law of seepage through the change of stress, that is, the coupling analysis of stress and seepage.Here, we take the length of the working face 220 m as the basis. In the mining geological conditions, we know that there is an aquifer above the bedrock, which will inevitably generate pore water pressure on the surrounding rock formations. The existence of pore water pressure will affect the surrounding rocks. The nature of the pore water will change, and the change of pore water pressure will affect the stress and strain of the rock formation, which will redistribute the stress of the rock formation. In the mining of the working face, the upper and lower ends of the working face and the cut-off are most affected. At this time, the permeability of the rock formation will increase accordingly. As the work surface advances, the stress of the overlying strata of the stope will be redistributed, and the overlying strata will be closed from the original tensile cracks and the permeability at this time The rate gradually decreases until the fissure is closed, the pore water pressure is zero, the working face continues to advance, and the pore water pressure may become negative, because the water-conducting fissure zone may spread to the overlying aquifer, affecting the water level of the aquifer.According to the simulation, the distribution of pore water pressure reflects the movement of the overlying strata in the process of advancing the working face. The figure shows that the stress distribution at both ends of the working face is concentrated, and a falling funnel appears above the goaf (Figure16).Figure 16
Cloud diagram of overburden pore water pressure at 250 m advance of working face advances to 100 m, 150 m, 200 m, and 250 m.
(a)
100 m advance of working face(b)
150 m advance of working face(c)
200 m advance of working face(d)
250 m advance of working faceAccording to the numerical simulation results, the change rule of overburden stress under the action of overburden pore water pressure can be obtained, which is expressed by the trend in Figure17.Figure 17
Comparison of permeability coefficient variation trend of different advancing distances at different working faces.
(a)(b)(c)Through the above analysis, the movement fracture of the overlying strata in coal mining is the result of the joint action of stress and seepage. The change of stress will change the flow law of fluid accordingly.


## 4. Conclusion

(1)
Summing up the relationship between stress and permeability coefficient in the earlier part, these relations are simple and convenient to apply(2)
By using computer numerical simulation technology, the influence of stress on the seepage of fractured rock mass is analyzed. The plastic separation diagram shows that the main areas of failure are above and on the left and right sides of the model. Under the influence of mining damage, the overburden failure is symmetrically distributed. By comprehensive comparison of the calculation models, different working face lengths lead to different distributions of plastic zone and development height of water-conducting fracture zone. Through the comparative analysis of cloud image and trend chart, it is found that the pore water pressure acts on the whole overlying strata and affects the development of fractures, and the change of fracture network also affects the subsidence trend of overlying aquifer, both of which influence each other(3)
The practice shows that the computer numerical simulation analysis of the coupling problem of stress and seepage is intuitive and effective. The application of clear and intuitive diagrams instead of complex theoretical formulas simplifies a lot of work(4)
The above analysis shows that the mining impact of the coal mining area in the Yellow River Basin has led to the destruction of the overlying rock aquifer and the leakage of surface water, affecting the fragile ecological environment of the surface, leading to lack of water for surface vegetation and soil desertification, so the mining area is underway. In ecological restoration, this method can be used to conduct in-depth research on the protection and restoration of aquifers, so as to provide a theoretical basis for the ecological protection and high-quality development of the Yellow River Basin



---

*Source: 4242337-2021-07-09.xml*