# Deflection Mechanism and Safety Analysis of Coal Mine Shaft in Deep Soil Strata

Authors: Jihuan Han; Jiuqun Zou; Chenchen Hu; Weihao Yang
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

The main shaft and auxiliary shaft in the Guotun Coal Mine underwent large deflections, with deflection values of 359 mm and 322 mm, respectively. These two deflections represent the first occurrence of such large vertical shaft deviations in the soil strata in China. The deflection problem has seriously affected the hoisting safety and lining safety and has become a serious impediment to the sustainable production of mines. Therefore, the deflection mechanism must be determined. For this purpose, based on mining subsidence theory, the spatial probability integral method and a more accurate time function were used to establish a model, called 3D dynamic prediction model, for predicting the shaft movement. The formulas for calculating the lining stress caused by coal mining were based on established models. With measured shaft deflection data, the prediction parameters for deep soil strata were calculated on the basis of an inversion analysis. A comparative analysis of measured and calculated deflection values revealed that the reason for shaft deflection in Guotun Coal Mine is the insufficient size of the protection coal pillar (PCP); namely, the design parameters of the PCP in current codes are not applicable to the deep soil strata. As a result, under the asymmetric mining conditions, mining causes the shaft to deflect without damage and under the symmetric mining conditions, mining causes the lining to fracture. The results have an extremely important significance for the prevention and control of shaft deflection, for the rational design of PCP, and for the sustainability of mine production.


---

## Body

## 1. Introduction

Since 2002, China has built 71 shafts with soil thickness over 400 m, and 86% of these shafts have a soil thickness greater than half of their depths (Figure1). The shaft is regarded as the throat of a mine, and its safety concerns the survival and sustainability of the entire mine. Since underground mining activities may have a certain effect on the shaft, the protection coal pillar (PCP) must be reserved. Unexpectedly, the shafts located in Guotun Coal Mine have experienced serious deflection disasters (Table 1 [1]). These deflections have caused the cage girders to deform and fall off. Considering that shaft deflection first occurred in China and is rare in the history of mining worldwide, the deflection mechanism must be clarified promptly. As is known to all, the shaft PCPs in European coal mines [2], as well as in China [3, 4], were mainly designed by considering the boundary angle or the movement angle. However, the angle values are similar under different mining conditions (Table 2). Combined with the actual situation of Guotun Coal Mine, the authors’ preliminary analysis indicated that the main reason for shaft deflection might be the insufficient size of PCP; that is, coal mining caused different displacements of the shaft at different depths. Based on the above analysis, mining subsidence theory was used to clarify the shaft deflection mechanism in deep soil strata. To study the problem of surface movement caused by mining, Litwiniszyn [5] proposed the stochastic medium theory. Later, Liu and Liao [6] further developed this theory into a probability integral method to predict mining subsidence. By combining the probability integral method and differential interferometric synthetic-aperture radar (SAR) technique, Fan et al. [7] established a model for extracting the large deformation mining subsidence; Diao et al. [8] proposed a new monitoring method that could retrieve highly accurate 3D displacement of mining subsidence; and Yang et al. [9] obtained an approach for cost-effective and accurate prediction of 3D mining-induced displacement under different extraction conditions. Using the probability integral method, Li et al. [10] calculated and simulated the influence of pillar mining subsidence on shaft safety. The key to the probability integral method is the solution of prediction parameters. For this purpose, Guo et al. [11] obtained more precise results by applying the theory of artificial neural network to this method, Zhang et al. [12] provided a highly effective method based on least squares support vector machine theory, and Han et al. [13] calculated the parameters that are suitable for deep soil strata by a back analysis. Regarding research on measuring and monitoring of subsidence, Zheng et al. [14] carried out monitoring and analysis of mining 3D deformation by multiplatform SAR images with the probability integral method. Considering that the above methods can predict only the final deformation, Knothe [15] built a time function model that can predict the dynamic surface subsidence based on the Mitscherlich growth law, which has widely been used and further improved [16–23]. Additionally, time functions were established by Liu et al. [24] based on the Harris model and by Xu and Li [25] based on the logistic growth model. However, due to many parameters, their models are relatively complex to solve.Figure 1
The proportion of soil thickness to shaft depth.Table 1
Shaft deflection values in Guotun Coal Mine.
Date (year. month)
Main shaft
Auxiliary shaft
Deflection values
mm
/
azimuth angle
°
2015.7
349
/
95
322
/
77
2017.8
359
/
79
318
/
75Table 2
Some soil movement angles in China.
Mine name
Coal thickness (m)
Mining depth (m)
Angles (°)
Jiaohe
1.0∼1.7
35∼110
45
Xuzhou
0.5∼5.59
90∼140
36∼45
Shuangyashan
0.8∼2.1
30∼220
45
Pingyuan
1.15∼3.54
100∼330
55
Huainan
1.8∼4.2
<180
40∼45
Weizhou
1.8∼3.5
<310
45∼55
Fuxin
1.5∼2.4
<400
40∼50
Fushun
20∼50
<540
45
Kailuan
0.9∼3.4
<600
35∼45
Zaozhuang
1.0∼1.7
<600
45
Jining
1.9∼4.3
>600
40∼45
Juye
0.9∼7.9
>600
40∼45Based on the Knothe time function, a more accurate time function was proposed. Combined with the probability integral method, a 3D dynamic prediction model was established to reveal the shaft deflection mechanism and to analyse the lining safety in deep soil strata. By comparing and analysing the theoretical and measured horizontal displacement values, the feasibility of the model was proved and the shaft deflection mechanism was clarified. The results in this paper are of great significance to the sustainability of mine production.


## 2. Models and Methods

### 2.1. Guotun Coal Mine

Guotun Coal Mine is located in Juye Coalfield in Shandong Province, with an annual production capacity of 2.4 million tons. All the shafts are located at the industrial square. The square is arranged in the coal-free zone. As of July 2018, the entire first mining area and part of the fourth area were already mined (Figure2). The start times of the first and fourth mining areas were January 2010 and October 2016, respectively. Based on the comprehensive analysis of the distribution, sizes, and mining times of the areas, it was concluded that the shafts in Guotun Coal Mine were mainly affected by the mining activities of the first mining area. The first mining area is approximately 6300 m long in the north and south directions and 3700 m wide in the east and west directions, with an average mining depth of 770 m and average coal thickness of 3.97 m. The coal seam is nearly horizontal, with a general dip angle of 5∼8°. According to the criterion [3], the square PCPs were designed by the vertical section method, with the soil and rock movement angles of 45° and 70°, respectively. The distance from the square boundary to the pillar boundary is approximately 718 m. The main shaft is 140 m from the square boundary in the west.Figure 2
The layout of mined working faces.


### 2.2. Time Function Solution

#### 2.2.1. Function Form

Assuming that the subsidence velocity at a certain moment is directly proportional to the difference between the dynamic subsidence value at this moment and the maximum subsidence value, Knothe [15] obtained a differential function as follows:(1)dwstdt=cwsm−wst,where wst is the instantaneous subsidence value of a point on the surface at time t, wsm is the maximum subsidence value of the point, and c is a time coefficient related to lithology.With the initial conditionws0=0, equation (1) can be derived as(2)wst=wsm1−e−ct.The Knothe time function can be expressed asTKnothet=1−e−ct. With the characteristics of the velocity and acceleration of the surface subsidence, the subsidence process should be a roughly S-shaped curve with time. Then, the time function can be improved by hypothesis to(3)Tt=1−1+λt−ct,where λ is an analogous coefficient to c.With proper values ofc and λ, Figure 3 shows that equation (3) is more accurate than the Knothe function in predicting subsidence.Figure 3
Influence of parameters on time function.


#### 2.2.2. Parameters Solution

The 1308 working face was first mined in July 2011, and the partial observation points arranged along the strike are shown in Figure4. The first survey was on July 29, 2011, and the second was on August 15, 2011: the maximum surface subsidence for two times was 6 mm and 340 mm (i.e., point Z16), respectively. To facilitate easy calculation, assume λ = 1. When taking the first survey time as the initial time, based on the measured subsidence values at points Z1 and Z2, the fitting curves and the corresponding parameter values are shown in Figure 5. The results showed that the values of c were feasible and it was reasonable to predict the subsidence of a point with λ = 1. If it is assumed that λ is a fixed value for other working faces, then c should be calculated according to the actual mining conditions.Figure 4
Locations of partial observation points in the 1308 working face.Figure 5
Subsidence curves at pointsZ1 and Z2.SupposeH is the mining depth, β is the main influence angle, r0 is the main influence radius at the surface, s0 is the displacement of the inflection point, v is the average mining rate, and lc is the goaf critical size, then r0 = H/tanβ, and lc = 2r0 + 2s0 (Figure 6). When the maximum subsidence is equal to 0.98wsm, the opening has reached its critical size [26], and the mining time t=lc/v. By combining equations (2) and (3), c can be derived as(4)c=−ln0.02tln1+t=−ln0.022vr0+s0ln1+2r0+2s0/v.Figure 6
Significance of subsidence parameters.


### 2.3. Model Solution for Stratal Movement

With coincident horizontal projection, a surface coordinate system and a mining coordinate system were established separately for each working face. When the parameters of thei-th working face are denoted by subscript i, the coordinates of mining unit B and point A are (ξi, ηi, Hi) and (xi, yi, zi), respectively (Figure 7). Point A′ and point B′ are the horizontal projections of point A and unit B at the surface. When the strike length of the i-th working face is li and the tendency length is si, based on the probability integral method and equation (3), the subsidence value of point A caused by mining the i-th working face at a certain moment is given as(5)wAixi,yi,zi,ti=w0i1−1+ti−citi×∫0li∫0siw—Bxi,yi,zidηidξi,w—Bxi,yi,zi=1rzi2exp−πxi−ξi2+yi−ηi2rzi2,where ti is the time of the i-th working face from the initial mining to a certain time, w0i is the maximum subsidence value, and rzi is the main influence radius at depth z, which can be written as(6)w0i=miqicosαi,rzi=r0iHi−ziHia=HitanβHi−ziHia,where mi and αi are the thickness and dip angle of the i-th working face, respectively, qi is the coefficient of stratal subsidence, and a is a constant.Figure 7
The spatial coordinate system. 1: surface; 2: working face.According to the superposition principle, the subsidence of pointA caused by mining n working faces at a certain moment can be given as(7)wA=∑i=1nwAixi,yi,zi,ti=∑i=1nw0i1−1+ti−citi×∫0li∫0si1rzi2exp−πxi−ξi2+yi−ηi2rzi2dηidξi.Equation (7) can be integrated as(8)wA=∑i=1nw0i1−1+ti−citi4erfkziyi−erfkziyi−si×erfkzixi−erfkzixi−li,where erf is the probability integral function, which can be calculated by(9)erfkzixi=2∫0kzixie−u2duπ,kzi=πrzi.As shown in Figure7, φi is the rotation angle from the xi-axis to the yi-axis in the positive direction and then to the specified direction at point A′. Assume χA is the inclination of point A in the φi direction, and then(10)χA=∑i=1nχAixi,yi,zi,ti,φi=∑i=1ncosφi∂wAixi,yi,zi,ti∂xi+sinφi∂wAixi,yi,zi,ti∂yi.According to reference [27], the horizontal displacement uA can be expressed as(11)uA=∑i=1nbrziχAixi,yi,zi,ti,φi,where b is the horizontal movement factor.Equation (11), a 3D dynamic prediction function, describes the horizontal movement of the strata caused by mining multiface. The detailed expression can be seen from equation (12), and the prediction parameters can be inversed by measured data. For those indeterminate parameters, the empirical values can be referred to [3].(12)uA=∑i=1nbw0i1−1+ti−citi2f1ixi,zig2iyi,zicosφi+f2iyi,zig1ixi,zisinφi,where(13)f1ixi,zi=exp−kzi2xi2−exp−ki2xi−li2,f2iyi,zi=exp−kzi2yi2−exp−ki2yi−si2,g1ixi,zi=erfkzixi−erfkzixi−li,g2iyi,zi=erfkziyi−erfkziyi−si.


### 2.4. Models Solution for Shaft Movement

Suppose the shaft depth isH0, the initial coordinates of the point at the position of shaft axis are (x0i, y0i, z0i) in the i-th surface coordinate system. Substituting the coordinates into equations (8) and (12), the vertical and horizontal displacements can be calculated. The maximum value of z0i should be smaller than the values of H0 and Hi. Considering that the shaft has a large slenderness ratio and a small lateral bending resistance, the displacement of the shaft is approximately equal to that of the soil at the same position. Taking the centre of the shaft head as the origin, then z = z0i must be true for any point on the shaft axis, and the coordinate system was established to describe the shaft horizontal bending (Figure 8). Using wz and u(z) to represent the horizontal and vertical displacements of those points on the shaft axis at different depths, equations (8) and (12) can be rewritten as(14)wz=∑i=1nw0i1−1+3ti−citi4g1ix0i,zg2iy0i,z,(15)uz=∑i=1nbw0i1−1+3ti−citi2f1ix0i,zg2iy0i,zcosφi+f2iy0i,zg1ix0i,zsinφi.Figure 8
The horizontal bending coordinate system.Under the mining of multiple faces, equations (14) and (15) are the 3D dynamic prediction models for the shaft movement in the vertical and horizontal directions.


### 2.5. Analysis of Shaft Lining Stress

When the coal was asymmetrically mined relative to the shaft, the strata produced vertical and horizontal displacements due to mining effects, wherein the vertical displacement produced downward additional stress on shaft lining, while the horizontal movement produced bending stress. Assuming that the lining is a homogeneous and continuous elastic material, based on the principle of superposition, the vertical and horizontal lining stresses can be analysed separately. Since the most likely occurrence position of lining rupture was at the inner edge [28], the vertical stress at this position was first analysed. Obviously, the vertical lining stress σz includes the self-weight stress σ0z, vertical additional stress σc(z, t), and bending stress σf(z, t). By assuming that the compressive stress is positive, σz can be expressed as(16)σz=σ0z+σcz,t±σfz,t.


#### 2.5.1. Self-Weight Stress

The self-weight stressσ0z can be determined as follows:(17)σ0z=γcz,where γc is the unit weight of the lining.


#### 2.5.2. Vertical Additional Stress

First, the vertical compressive strainεc(z, t) of the lining can be calculated by equation (14), which can be derived as(18)εcz=∑i=1nw0i1−1+ti−citi4g2iy0i,z∂g1ix0i,z∂z+g1ix0i,z∂g2iy0i,z∂z.Then, the vertical additional stressσc(z, t) can be obtained with equation (18) as(19)σcz,t=−Ecεcz,t,where Ec is the elastic modulus of the lining and “−” is due to the positive compressive stress.


#### 2.5.3. Bending Stress

When the shaft can be regarded as an elastic beam, equation (15) can be expressed as the flexural formula of the beam. Based on the theory of the elastic beam, the bending moment of the lining can be obtained as(20)Mz,t=−EcIzd2uz,tdz2,(21)Iz=πDw4−d464,where Iz is the moment of inertia and d and Dw are the net diameter of the shaft and the outer diameter of the outer lining, respectively.Therefore, the bending stressσf(z, t) at the inner edge of the lining can be obtained based on equations (20) and (21) as follows:(22)σfz,t=32dMz,tπDw4−d4.Compared to asymmetric mining, the strata produced double vertical displacement to the shaft under the symmetrical mining. The vertical lining stressσz′ only included the self-weight stress σ0z and the vertical additional stress σc′z,t, here σc′z,t=2σcz,t. Then, σz′ can be expressed as(23)σz′=σ0z+2σcz,t.


### 2.6. Analysis of Lining Safety

According to the theory of a thick cylindrical lining, the circumferential stressσθ and radial stress σr at the inner edge of the lining can be calculated as(24)σr=0,σθ=2pwDn2Dn2−d2,where pw is the hydrostatic pressure and Dn is the outer diameter of the inner shaft lining.According to Code [29], when equation (25) is satisfied, the inner edge of the lining is in an unsafe state under the influence of coal mining:(25)σz>ηfc+ρminfy′,σz′>ηfc+ρminfy′,where fc is the uniaxial compressive design strength of concrete, fy′ is the compressive design strength of reinforcement, η is the improvement coefficient of the concrete strength under multiaxial stress, η is related to σθ and σr, and in general, η = 1.2. Finally, ρmin is the minimum steel content.


## 2.1. Guotun Coal Mine

Guotun Coal Mine is located in Juye Coalfield in Shandong Province, with an annual production capacity of 2.4 million tons. All the shafts are located at the industrial square. The square is arranged in the coal-free zone. As of July 2018, the entire first mining area and part of the fourth area were already mined (Figure2). The start times of the first and fourth mining areas were January 2010 and October 2016, respectively. Based on the comprehensive analysis of the distribution, sizes, and mining times of the areas, it was concluded that the shafts in Guotun Coal Mine were mainly affected by the mining activities of the first mining area. The first mining area is approximately 6300 m long in the north and south directions and 3700 m wide in the east and west directions, with an average mining depth of 770 m and average coal thickness of 3.97 m. The coal seam is nearly horizontal, with a general dip angle of 5∼8°. According to the criterion [3], the square PCPs were designed by the vertical section method, with the soil and rock movement angles of 45° and 70°, respectively. The distance from the square boundary to the pillar boundary is approximately 718 m. The main shaft is 140 m from the square boundary in the west.Figure 2
The layout of mined working faces.


## 2.2. Time Function Solution

### 2.2.1. Function Form

Assuming that the subsidence velocity at a certain moment is directly proportional to the difference between the dynamic subsidence value at this moment and the maximum subsidence value, Knothe [15] obtained a differential function as follows:(1)dwstdt=cwsm−wst,where wst is the instantaneous subsidence value of a point on the surface at time t, wsm is the maximum subsidence value of the point, and c is a time coefficient related to lithology.With the initial conditionws0=0, equation (1) can be derived as(2)wst=wsm1−e−ct.The Knothe time function can be expressed asTKnothet=1−e−ct. With the characteristics of the velocity and acceleration of the surface subsidence, the subsidence process should be a roughly S-shaped curve with time. Then, the time function can be improved by hypothesis to(3)Tt=1−1+λt−ct,where λ is an analogous coefficient to c.With proper values ofc and λ, Figure 3 shows that equation (3) is more accurate than the Knothe function in predicting subsidence.Figure 3
Influence of parameters on time function.


### 2.2.2. Parameters Solution

The 1308 working face was first mined in July 2011, and the partial observation points arranged along the strike are shown in Figure4. The first survey was on July 29, 2011, and the second was on August 15, 2011: the maximum surface subsidence for two times was 6 mm and 340 mm (i.e., point Z16), respectively. To facilitate easy calculation, assume λ = 1. When taking the first survey time as the initial time, based on the measured subsidence values at points Z1 and Z2, the fitting curves and the corresponding parameter values are shown in Figure 5. The results showed that the values of c were feasible and it was reasonable to predict the subsidence of a point with λ = 1. If it is assumed that λ is a fixed value for other working faces, then c should be calculated according to the actual mining conditions.Figure 4
Locations of partial observation points in the 1308 working face.Figure 5
Subsidence curves at pointsZ1 and Z2.SupposeH is the mining depth, β is the main influence angle, r0 is the main influence radius at the surface, s0 is the displacement of the inflection point, v is the average mining rate, and lc is the goaf critical size, then r0 = H/tanβ, and lc = 2r0 + 2s0 (Figure 6). When the maximum subsidence is equal to 0.98wsm, the opening has reached its critical size [26], and the mining time t=lc/v. By combining equations (2) and (3), c can be derived as(4)c=−ln0.02tln1+t=−ln0.022vr0+s0ln1+2r0+2s0/v.Figure 6
Significance of subsidence parameters.


## 2.2.1. Function Form

Assuming that the subsidence velocity at a certain moment is directly proportional to the difference between the dynamic subsidence value at this moment and the maximum subsidence value, Knothe [15] obtained a differential function as follows:(1)dwstdt=cwsm−wst,where wst is the instantaneous subsidence value of a point on the surface at time t, wsm is the maximum subsidence value of the point, and c is a time coefficient related to lithology.With the initial conditionws0=0, equation (1) can be derived as(2)wst=wsm1−e−ct.The Knothe time function can be expressed asTKnothet=1−e−ct. With the characteristics of the velocity and acceleration of the surface subsidence, the subsidence process should be a roughly S-shaped curve with time. Then, the time function can be improved by hypothesis to(3)Tt=1−1+λt−ct,where λ is an analogous coefficient to c.With proper values ofc and λ, Figure 3 shows that equation (3) is more accurate than the Knothe function in predicting subsidence.Figure 3
Influence of parameters on time function.


## 2.2.2. Parameters Solution

The 1308 working face was first mined in July 2011, and the partial observation points arranged along the strike are shown in Figure4. The first survey was on July 29, 2011, and the second was on August 15, 2011: the maximum surface subsidence for two times was 6 mm and 340 mm (i.e., point Z16), respectively. To facilitate easy calculation, assume λ = 1. When taking the first survey time as the initial time, based on the measured subsidence values at points Z1 and Z2, the fitting curves and the corresponding parameter values are shown in Figure 5. The results showed that the values of c were feasible and it was reasonable to predict the subsidence of a point with λ = 1. If it is assumed that λ is a fixed value for other working faces, then c should be calculated according to the actual mining conditions.Figure 4
Locations of partial observation points in the 1308 working face.Figure 5
Subsidence curves at pointsZ1 and Z2.SupposeH is the mining depth, β is the main influence angle, r0 is the main influence radius at the surface, s0 is the displacement of the inflection point, v is the average mining rate, and lc is the goaf critical size, then r0 = H/tanβ, and lc = 2r0 + 2s0 (Figure 6). When the maximum subsidence is equal to 0.98wsm, the opening has reached its critical size [26], and the mining time t=lc/v. By combining equations (2) and (3), c can be derived as(4)c=−ln0.02tln1+t=−ln0.022vr0+s0ln1+2r0+2s0/v.Figure 6
Significance of subsidence parameters.


## 2.3. Model Solution for Stratal Movement

With coincident horizontal projection, a surface coordinate system and a mining coordinate system were established separately for each working face. When the parameters of thei-th working face are denoted by subscript i, the coordinates of mining unit B and point A are (ξi, ηi, Hi) and (xi, yi, zi), respectively (Figure 7). Point A′ and point B′ are the horizontal projections of point A and unit B at the surface. When the strike length of the i-th working face is li and the tendency length is si, based on the probability integral method and equation (3), the subsidence value of point A caused by mining the i-th working face at a certain moment is given as(5)wAixi,yi,zi,ti=w0i1−1+ti−citi×∫0li∫0siw—Bxi,yi,zidηidξi,w—Bxi,yi,zi=1rzi2exp−πxi−ξi2+yi−ηi2rzi2,where ti is the time of the i-th working face from the initial mining to a certain time, w0i is the maximum subsidence value, and rzi is the main influence radius at depth z, which can be written as(6)w0i=miqicosαi,rzi=r0iHi−ziHia=HitanβHi−ziHia,where mi and αi are the thickness and dip angle of the i-th working face, respectively, qi is the coefficient of stratal subsidence, and a is a constant.Figure 7
The spatial coordinate system. 1: surface; 2: working face.According to the superposition principle, the subsidence of pointA caused by mining n working faces at a certain moment can be given as(7)wA=∑i=1nwAixi,yi,zi,ti=∑i=1nw0i1−1+ti−citi×∫0li∫0si1rzi2exp−πxi−ξi2+yi−ηi2rzi2dηidξi.Equation (7) can be integrated as(8)wA=∑i=1nw0i1−1+ti−citi4erfkziyi−erfkziyi−si×erfkzixi−erfkzixi−li,where erf is the probability integral function, which can be calculated by(9)erfkzixi=2∫0kzixie−u2duπ,kzi=πrzi.As shown in Figure7, φi is the rotation angle from the xi-axis to the yi-axis in the positive direction and then to the specified direction at point A′. Assume χA is the inclination of point A in the φi direction, and then(10)χA=∑i=1nχAixi,yi,zi,ti,φi=∑i=1ncosφi∂wAixi,yi,zi,ti∂xi+sinφi∂wAixi,yi,zi,ti∂yi.According to reference [27], the horizontal displacement uA can be expressed as(11)uA=∑i=1nbrziχAixi,yi,zi,ti,φi,where b is the horizontal movement factor.Equation (11), a 3D dynamic prediction function, describes the horizontal movement of the strata caused by mining multiface. The detailed expression can be seen from equation (12), and the prediction parameters can be inversed by measured data. For those indeterminate parameters, the empirical values can be referred to [3].(12)uA=∑i=1nbw0i1−1+ti−citi2f1ixi,zig2iyi,zicosφi+f2iyi,zig1ixi,zisinφi,where(13)f1ixi,zi=exp−kzi2xi2−exp−ki2xi−li2,f2iyi,zi=exp−kzi2yi2−exp−ki2yi−si2,g1ixi,zi=erfkzixi−erfkzixi−li,g2iyi,zi=erfkziyi−erfkziyi−si.


## 2.4. Models Solution for Shaft Movement

Suppose the shaft depth isH0, the initial coordinates of the point at the position of shaft axis are (x0i, y0i, z0i) in the i-th surface coordinate system. Substituting the coordinates into equations (8) and (12), the vertical and horizontal displacements can be calculated. The maximum value of z0i should be smaller than the values of H0 and Hi. Considering that the shaft has a large slenderness ratio and a small lateral bending resistance, the displacement of the shaft is approximately equal to that of the soil at the same position. Taking the centre of the shaft head as the origin, then z = z0i must be true for any point on the shaft axis, and the coordinate system was established to describe the shaft horizontal bending (Figure 8). Using wz and u(z) to represent the horizontal and vertical displacements of those points on the shaft axis at different depths, equations (8) and (12) can be rewritten as(14)wz=∑i=1nw0i1−1+3ti−citi4g1ix0i,zg2iy0i,z,(15)uz=∑i=1nbw0i1−1+3ti−citi2f1ix0i,zg2iy0i,zcosφi+f2iy0i,zg1ix0i,zsinφi.Figure 8
The horizontal bending coordinate system.Under the mining of multiple faces, equations (14) and (15) are the 3D dynamic prediction models for the shaft movement in the vertical and horizontal directions.


## 2.5. Analysis of Shaft Lining Stress

When the coal was asymmetrically mined relative to the shaft, the strata produced vertical and horizontal displacements due to mining effects, wherein the vertical displacement produced downward additional stress on shaft lining, while the horizontal movement produced bending stress. Assuming that the lining is a homogeneous and continuous elastic material, based on the principle of superposition, the vertical and horizontal lining stresses can be analysed separately. Since the most likely occurrence position of lining rupture was at the inner edge [28], the vertical stress at this position was first analysed. Obviously, the vertical lining stress σz includes the self-weight stress σ0z, vertical additional stress σc(z, t), and bending stress σf(z, t). By assuming that the compressive stress is positive, σz can be expressed as(16)σz=σ0z+σcz,t±σfz,t.


### 2.5.1. Self-Weight Stress

The self-weight stressσ0z can be determined as follows:(17)σ0z=γcz,where γc is the unit weight of the lining.


### 2.5.2. Vertical Additional Stress

First, the vertical compressive strainεc(z, t) of the lining can be calculated by equation (14), which can be derived as(18)εcz=∑i=1nw0i1−1+ti−citi4g2iy0i,z∂g1ix0i,z∂z+g1ix0i,z∂g2iy0i,z∂z.Then, the vertical additional stressσc(z, t) can be obtained with equation (18) as(19)σcz,t=−Ecεcz,t,where Ec is the elastic modulus of the lining and “−” is due to the positive compressive stress.


### 2.5.3. Bending Stress

When the shaft can be regarded as an elastic beam, equation (15) can be expressed as the flexural formula of the beam. Based on the theory of the elastic beam, the bending moment of the lining can be obtained as(20)Mz,t=−EcIzd2uz,tdz2,(21)Iz=πDw4−d464,where Iz is the moment of inertia and d and Dw are the net diameter of the shaft and the outer diameter of the outer lining, respectively.Therefore, the bending stressσf(z, t) at the inner edge of the lining can be obtained based on equations (20) and (21) as follows:(22)σfz,t=32dMz,tπDw4−d4.Compared to asymmetric mining, the strata produced double vertical displacement to the shaft under the symmetrical mining. The vertical lining stressσz′ only included the self-weight stress σ0z and the vertical additional stress σc′z,t, here σc′z,t=2σcz,t. Then, σz′ can be expressed as(23)σz′=σ0z+2σcz,t.


## 2.5.1. Self-Weight Stress

The self-weight stressσ0z can be determined as follows:(17)σ0z=γcz,where γc is the unit weight of the lining.


## 2.5.2. Vertical Additional Stress

First, the vertical compressive strainεc(z, t) of the lining can be calculated by equation (14), which can be derived as(18)εcz=∑i=1nw0i1−1+ti−citi4g2iy0i,z∂g1ix0i,z∂z+g1ix0i,z∂g2iy0i,z∂z.Then, the vertical additional stressσc(z, t) can be obtained with equation (18) as(19)σcz,t=−Ecεcz,t,where Ec is the elastic modulus of the lining and “−” is due to the positive compressive stress.


## 2.5.3. Bending Stress

When the shaft can be regarded as an elastic beam, equation (15) can be expressed as the flexural formula of the beam. Based on the theory of the elastic beam, the bending moment of the lining can be obtained as(20)Mz,t=−EcIzd2uz,tdz2,(21)Iz=πDw4−d464,where Iz is the moment of inertia and d and Dw are the net diameter of the shaft and the outer diameter of the outer lining, respectively.Therefore, the bending stressσf(z, t) at the inner edge of the lining can be obtained based on equations (20) and (21) as follows:(22)σfz,t=32dMz,tπDw4−d4.Compared to asymmetric mining, the strata produced double vertical displacement to the shaft under the symmetrical mining. The vertical lining stressσz′ only included the self-weight stress σ0z and the vertical additional stress σc′z,t, here σc′z,t=2σcz,t. Then, σz′ can be expressed as(23)σz′=σ0z+2σcz,t.


## 2.6. Analysis of Lining Safety

According to the theory of a thick cylindrical lining, the circumferential stressσθ and radial stress σr at the inner edge of the lining can be calculated as(24)σr=0,σθ=2pwDn2Dn2−d2,where pw is the hydrostatic pressure and Dn is the outer diameter of the inner shaft lining.According to Code [29], when equation (25) is satisfied, the inner edge of the lining is in an unsafe state under the influence of coal mining:(25)σz>ηfc+ρminfy′,σz′>ηfc+ρminfy′,where fc is the uniaxial compressive design strength of concrete, fy′ is the compressive design strength of reinforcement, η is the improvement coefficient of the concrete strength under multiaxial stress, η is related to σθ and σr, and in general, η = 1.2. Finally, ρmin is the minimum steel content.


## 3. Analysis of Deflection Mechanism and Lining Safety

### 3.1. Deflection Mechanism Analysis

#### 3.1.1. Prediction Parameters Solution

The prediction parameters in Guotun Coal Mine were obtained as follows:s0i = 0 [3] and q = 1.0 [30], and other parameters were taken based on an inversion analysis. The horizontal displacements of the main shaft were expressed by um. An inversion analysis was carried out by the measured values of um in June 2015 (curve I in Figure 9), and the inversion results were basically consistent with the measured values and were suitable for equation (15). Using the inversion parameters, the values of um in 2017 were calculated as seen from curve II in Figure 9. Compared with the measured values, their average difference was only 2.1 mm and smaller than the measurement error (8 mm). With high precision of anastomosis, the horizontal shaft displacement at different depths could be inverted by the dynamic prediction model. According to the calculation results, the prediction parameters obtained by inversion were a = 0.7, b = 0.66, and tanβ = 0.50. The calculation parameters of the working faces and shaft coordinates are shown in Table 3.Figure 9
The horizontal displacement of the main shaft.Table 3
Corresponding parameter values in Guotun Coal Mine.
i
Face name
Main shaft
Auxiliary shaft
Start date (year. month)
End data (year. month)
H
i (m)
m
i (m)
φ
i (°)
s
i (m)
l
i (m)
x
i (m)
y
i (m)
x
i (m)
y
i (m)
1
1301
2525
692
2445
737
2010.11
2013.01
840
3.17
335
227
1768
2
1302
810
950
855
1030
2010.01
2010.07
780
2.80
178
200
490
3
1303
2535
937
2455
978
2012.06
2015.06
845
2.99
335
230
1675
2015.08
1695
4
1304
785
1378
826
1460
2010.08
2011.03
770
3.20
202
150
870
5
1305
2531
1194
2451
1236
2013.09
2015.06
820
3.41
335
245
985
2016.08
1250
6
1306
942
1913
1016
1968
2015.07
2015.12
745
2.80
208
110
460
7
1307
2836
1445
2816
1490
2016.04
2017.08
770
3.30
335
240
2100
8
1308
1280
1385
1330
1462
2011.07
2012.03
740
3.08
168
157
888
9
1309
2868
1630
2039
1708
2017.08
2018.07
700
3.30
335
170
1200
10
1310
796
2065
836
2148
2012.12
2013.08
730
3.20
198
190
790
11
1311
2060
349
2136
297
2017.04
2017.07
785
3.10
105
136
266
12
1312-1
655
2313
677
2402
2014.03
2014.06
730
3.00
208
130
640
13
1312-2
934
2137
981
2216
2014.07
2015.05
735
3.00
188
130
1340
14
1314
957
2499
311
2731
2017.12
2018.07
710
2.60
189
199
690
15
1315
1090
1605
1177
1635
2016.07
2016.09
740
3.10
169
110
360


#### 3.1.2. Prediction Effect Verification

To further verify the feasibility of the prediction parameters, the auxiliary shaft was taken as a case. The horizontal displacements of the auxiliary shaft were expressed byua. The theoretical values and the measured values in 2015 (curve I) and 2017 (curve II) are shown in Figure 10. The errors in two times were only approximately 7% and 2% of the measured values at the shaft head, respectively, which indicated that the calculation results were still highly reliable considering the measurement error and the difference between the main shaft and auxiliary shaft.Figure 10
The horizontal displacement of the auxiliary shaft.


#### 3.1.3. Analysis of Results

The prediction model can accurately predict the shaft deflection in deep soil strata and calculate the horizontal displacement magnitudes at different depths under conditions of multiface mining. The calculated results agreed with the measured values. Based on the prediction parameters of inversion, it could be calculated that whenH > 700 m, r0 = H/tan β > 1400 m. Therefore, the reason for the shaft deflection in Guotun Coal Mine was the insufficient size of PCP.


### 3.2. Lining Safety Analysis

Taking the main shaft as an example, the vertical theoretical displacement was calculated using the 3D dynamic prediction model (Figure11). The calculation results showed that the maximum vertical displacement at two times was 37.3% and 35.4% of the maximum horizontal displacement, indicating that the mining effect was greater on the horizontal displacement of the shaft than on the vertical displacement in deep soil strata.Figure 11
Vertical theoretical displacement of main shaft lining caused by coal mining.The azimuth angles of the main shaft at different depths in 2015 were calculated by the measured and theoretical displacement values, respectively. The azimuth angles did not change with depth in the soil section, and the two calculation results were basically consistent (Figure12). When using the 3D dynamic prediction model to calculate the force of the shaft in the soil section, it was reasonable to consider the main shaft as a plane bending deformation.Figure 12
The azimuth angle of main shaft deflection.To comprehensively analyse the influence of mining on lining safety, the vertical stress was analysed under symmetrical and asymmetric mining conditions in Guotun Coal Mine. According to the actual parameters of the lining structure of the main shaft in the soil section (Table4), the vertical stress in 2015 and 2017 was calculated. The main shaft lining adopted an HRB335 steel bar, for which fy′=300MPa. When ρmin = 0.2%, the stress values at different times are shown in Figure 13.Table 4
Actual parameters of shaft lining in the soil section.
z (m)
Strength grade
E
c (×104 MPa)
f
c (MPa)
D
n (mm)
D
w (mm)
0∼100
C30
3.00
14.3
6100
7200
100∼210
C40
3.25
19.1
210∼247
C40
3.25
19.1
6400
7800
247∼300
C50
3.45
23.1
300∼340
C55
3.55
25.3
340∼370
C60
3.60
27.5
370∼400
C50
3.45
23.1
7000
9000
400∼446
C55
3.55
25.3
446∼498
C60
3.60
27.5
498∼545
C65
3.65
29.7
9300
545∼595
C70
3.70
31.8Figure 13
Safety analysis of main shaft lining.Figure13 shows that, under the conditions of asymmetric mining, the vertical stress values are less than the yield strength values, and the main shaft remains safe. Shaft deflection mainly affects hoisting safety and has a smaller effect on lining safety. However, if the coal is mined symmetrically in Guotun Coal Mine, the mining will cause the lining to fracture. In conclusion, the problems are mainly caused by unreasonable design size of PCP in deep soil strata.


## 3.1. Deflection Mechanism Analysis

### 3.1.1. Prediction Parameters Solution

The prediction parameters in Guotun Coal Mine were obtained as follows:s0i = 0 [3] and q = 1.0 [30], and other parameters were taken based on an inversion analysis. The horizontal displacements of the main shaft were expressed by um. An inversion analysis was carried out by the measured values of um in June 2015 (curve I in Figure 9), and the inversion results were basically consistent with the measured values and were suitable for equation (15). Using the inversion parameters, the values of um in 2017 were calculated as seen from curve II in Figure 9. Compared with the measured values, their average difference was only 2.1 mm and smaller than the measurement error (8 mm). With high precision of anastomosis, the horizontal shaft displacement at different depths could be inverted by the dynamic prediction model. According to the calculation results, the prediction parameters obtained by inversion were a = 0.7, b = 0.66, and tanβ = 0.50. The calculation parameters of the working faces and shaft coordinates are shown in Table 3.Figure 9
The horizontal displacement of the main shaft.Table 3
Corresponding parameter values in Guotun Coal Mine.
i
Face name
Main shaft
Auxiliary shaft
Start date (year. month)
End data (year. month)
H
i (m)
m
i (m)
φ
i (°)
s
i (m)
l
i (m)
x
i (m)
y
i (m)
x
i (m)
y
i (m)
1
1301
2525
692
2445
737
2010.11
2013.01
840
3.17
335
227
1768
2
1302
810
950
855
1030
2010.01
2010.07
780
2.80
178
200
490
3
1303
2535
937
2455
978
2012.06
2015.06
845
2.99
335
230
1675
2015.08
1695
4
1304
785
1378
826
1460
2010.08
2011.03
770
3.20
202
150
870
5
1305
2531
1194
2451
1236
2013.09
2015.06
820
3.41
335
245
985
2016.08
1250
6
1306
942
1913
1016
1968
2015.07
2015.12
745
2.80
208
110
460
7
1307
2836
1445
2816
1490
2016.04
2017.08
770
3.30
335
240
2100
8
1308
1280
1385
1330
1462
2011.07
2012.03
740
3.08
168
157
888
9
1309
2868
1630
2039
1708
2017.08
2018.07
700
3.30
335
170
1200
10
1310
796
2065
836
2148
2012.12
2013.08
730
3.20
198
190
790
11
1311
2060
349
2136
297
2017.04
2017.07
785
3.10
105
136
266
12
1312-1
655
2313
677
2402
2014.03
2014.06
730
3.00
208
130
640
13
1312-2
934
2137
981
2216
2014.07
2015.05
735
3.00
188
130
1340
14
1314
957
2499
311
2731
2017.12
2018.07
710
2.60
189
199
690
15
1315
1090
1605
1177
1635
2016.07
2016.09
740
3.10
169
110
360


### 3.1.2. Prediction Effect Verification

To further verify the feasibility of the prediction parameters, the auxiliary shaft was taken as a case. The horizontal displacements of the auxiliary shaft were expressed byua. The theoretical values and the measured values in 2015 (curve I) and 2017 (curve II) are shown in Figure 10. The errors in two times were only approximately 7% and 2% of the measured values at the shaft head, respectively, which indicated that the calculation results were still highly reliable considering the measurement error and the difference between the main shaft and auxiliary shaft.Figure 10
The horizontal displacement of the auxiliary shaft.


### 3.1.3. Analysis of Results

The prediction model can accurately predict the shaft deflection in deep soil strata and calculate the horizontal displacement magnitudes at different depths under conditions of multiface mining. The calculated results agreed with the measured values. Based on the prediction parameters of inversion, it could be calculated that whenH > 700 m, r0 = H/tan β > 1400 m. Therefore, the reason for the shaft deflection in Guotun Coal Mine was the insufficient size of PCP.


## 3.1.1. Prediction Parameters Solution

The prediction parameters in Guotun Coal Mine were obtained as follows:s0i = 0 [3] and q = 1.0 [30], and other parameters were taken based on an inversion analysis. The horizontal displacements of the main shaft were expressed by um. An inversion analysis was carried out by the measured values of um in June 2015 (curve I in Figure 9), and the inversion results were basically consistent with the measured values and were suitable for equation (15). Using the inversion parameters, the values of um in 2017 were calculated as seen from curve II in Figure 9. Compared with the measured values, their average difference was only 2.1 mm and smaller than the measurement error (8 mm). With high precision of anastomosis, the horizontal shaft displacement at different depths could be inverted by the dynamic prediction model. According to the calculation results, the prediction parameters obtained by inversion were a = 0.7, b = 0.66, and tanβ = 0.50. The calculation parameters of the working faces and shaft coordinates are shown in Table 3.Figure 9
The horizontal displacement of the main shaft.Table 3
Corresponding parameter values in Guotun Coal Mine.
i
Face name
Main shaft
Auxiliary shaft
Start date (year. month)
End data (year. month)
H
i (m)
m
i (m)
φ
i (°)
s
i (m)
l
i (m)
x
i (m)
y
i (m)
x
i (m)
y
i (m)
1
1301
2525
692
2445
737
2010.11
2013.01
840
3.17
335
227
1768
2
1302
810
950
855
1030
2010.01
2010.07
780
2.80
178
200
490
3
1303
2535
937
2455
978
2012.06
2015.06
845
2.99
335
230
1675
2015.08
1695
4
1304
785
1378
826
1460
2010.08
2011.03
770
3.20
202
150
870
5
1305
2531
1194
2451
1236
2013.09
2015.06
820
3.41
335
245
985
2016.08
1250
6
1306
942
1913
1016
1968
2015.07
2015.12
745
2.80
208
110
460
7
1307
2836
1445
2816
1490
2016.04
2017.08
770
3.30
335
240
2100
8
1308
1280
1385
1330
1462
2011.07
2012.03
740
3.08
168
157
888
9
1309
2868
1630
2039
1708
2017.08
2018.07
700
3.30
335
170
1200
10
1310
796
2065
836
2148
2012.12
2013.08
730
3.20
198
190
790
11
1311
2060
349
2136
297
2017.04
2017.07
785
3.10
105
136
266
12
1312-1
655
2313
677
2402
2014.03
2014.06
730
3.00
208
130
640
13
1312-2
934
2137
981
2216
2014.07
2015.05
735
3.00
188
130
1340
14
1314
957
2499
311
2731
2017.12
2018.07
710
2.60
189
199
690
15
1315
1090
1605
1177
1635
2016.07
2016.09
740
3.10
169
110
360


## 3.1.2. Prediction Effect Verification

To further verify the feasibility of the prediction parameters, the auxiliary shaft was taken as a case. The horizontal displacements of the auxiliary shaft were expressed byua. The theoretical values and the measured values in 2015 (curve I) and 2017 (curve II) are shown in Figure 10. The errors in two times were only approximately 7% and 2% of the measured values at the shaft head, respectively, which indicated that the calculation results were still highly reliable considering the measurement error and the difference between the main shaft and auxiliary shaft.Figure 10
The horizontal displacement of the auxiliary shaft.


## 3.1.3. Analysis of Results

The prediction model can accurately predict the shaft deflection in deep soil strata and calculate the horizontal displacement magnitudes at different depths under conditions of multiface mining. The calculated results agreed with the measured values. Based on the prediction parameters of inversion, it could be calculated that whenH > 700 m, r0 = H/tan β > 1400 m. Therefore, the reason for the shaft deflection in Guotun Coal Mine was the insufficient size of PCP.


## 3.2. Lining Safety Analysis

Taking the main shaft as an example, the vertical theoretical displacement was calculated using the 3D dynamic prediction model (Figure11). The calculation results showed that the maximum vertical displacement at two times was 37.3% and 35.4% of the maximum horizontal displacement, indicating that the mining effect was greater on the horizontal displacement of the shaft than on the vertical displacement in deep soil strata.Figure 11
Vertical theoretical displacement of main shaft lining caused by coal mining.The azimuth angles of the main shaft at different depths in 2015 were calculated by the measured and theoretical displacement values, respectively. The azimuth angles did not change with depth in the soil section, and the two calculation results were basically consistent (Figure12). When using the 3D dynamic prediction model to calculate the force of the shaft in the soil section, it was reasonable to consider the main shaft as a plane bending deformation.Figure 12
The azimuth angle of main shaft deflection.To comprehensively analyse the influence of mining on lining safety, the vertical stress was analysed under symmetrical and asymmetric mining conditions in Guotun Coal Mine. According to the actual parameters of the lining structure of the main shaft in the soil section (Table4), the vertical stress in 2015 and 2017 was calculated. The main shaft lining adopted an HRB335 steel bar, for which fy′=300MPa. When ρmin = 0.2%, the stress values at different times are shown in Figure 13.Table 4
Actual parameters of shaft lining in the soil section.
z (m)
Strength grade
E
c (×104 MPa)
f
c (MPa)
D
n (mm)
D
w (mm)
0∼100
C30
3.00
14.3
6100
7200
100∼210
C40
3.25
19.1
210∼247
C40
3.25
19.1
6400
7800
247∼300
C50
3.45
23.1
300∼340
C55
3.55
25.3
340∼370
C60
3.60
27.5
370∼400
C50
3.45
23.1
7000
9000
400∼446
C55
3.55
25.3
446∼498
C60
3.60
27.5
498∼545
C65
3.65
29.7
9300
545∼595
C70
3.70
31.8Figure 13
Safety analysis of main shaft lining.Figure13 shows that, under the conditions of asymmetric mining, the vertical stress values are less than the yield strength values, and the main shaft remains safe. Shaft deflection mainly affects hoisting safety and has a smaller effect on lining safety. However, if the coal is mined symmetrically in Guotun Coal Mine, the mining will cause the lining to fracture. In conclusion, the problems are mainly caused by unreasonable design size of PCP in deep soil strata.


## 4. Conclusions

Combining the probability integral method with a new time function, a 3D dynamic prediction model was established to describe the shaft movement. Based on the shaft deflection in Guotun Coal Mine, the prediction parameters suitable for deep soil strata were obtained. By comparing and analysing the theoretical and measured values, the following conclusions were obtained:(1)
The shaft deflection mechanism in deep soil strata was as follows: under the conditions of asymmetric mining, the mining caused the soil around the shaft to produce different displacement magnitudes at different depths, resulting in the vertical and horizontal displacements of the shaft within the range of mining influence. According to the current regulations, the shafts in Guotun Coal Mine are not located in the scope of mining influence, while they have experienced severe deflections. These deflections indicate that the stipulated design parameters for the shaft PCPs in current codes are not applicable to the deep soil strata.(2)
The new time function was more accurate in predicting subsidence than the Knothe function. Withλ = 1, the parameters of the prediction model were obtained by inverting the measured data. The prediction model can accurately predict the shaft movement, which provides an important reference for preventing shaft deflection and designing the shaft PCPs in deep soil strata.(3)
Considering the current situation in Guotun Coal Mine, the grouting technology should be adopted to correct the deflection. For similar newly built shafts, the symmetrical or filling mining is recommended according to the actual situation, and lining structures should adapt themselves to a certain bending deformation and be installed with contractible devices.



---

*Source: 9461742-2019-12-05.xml*