# Study on the Thin Plate Model with Elastic Foundation Boundary of Overlying Strata for Backfill Mining

Authors: Dongdong Chen; Xiaoyu Wu; Shengrong Xie; Yanding Sun; Qing Zhang; En Wang; Yaohui Sun; Long Wang; Hui Li; Zaisheng Jiang; Xiaowei Wu
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

In order to better study the movement principles of overlying strata during backfill mining, we established a thin plate model on an elastic foundation with elastic foundation boundary of the main roof. And by the finite difference method, the variation principles of the main roof’s principal moments and maximum subsidenceω0 with the elastic foundation coefficient k1 of the coal seam, the elastic foundation coefficient k2 of backfill body, the thickness h, Young’s modulus E, and Poisson’s ratio μ of main roof are calculated and studied. Using these calculations, we were able to determine that the main roof had three principal bending moment extreme points, including Mzz in backfill areas, Mc of the long side area, and Md of the short side area. The distance Lc of Mc advancing coal wall continuously increased with the increase in k2, while the principal moment of main roof’s middle area decreased with an increase in k2; when k2 became larger, the maximum principal moment in the midpoint of main roof transferred to the surrounding and the maximum principal moments was in four-corner area; Mc and Md decreased with an increase in k2, and Md was more sensitive to k2 than Mc; and Md decreased significantly with the increase in k2. Lc continuously decreased with the increase in k1, while Mc, Md and Mzz increased with the increase in k1 and the reduced amplitude of Mzz was the minimum. The effect of μ on principal bending moments and ω0 was very small; The growth rate of Mzz was the largest when E or h increased. Md, Mzz, and Lc remained unchanged when k1, k2, and Young’s modulus E of the main roof increased while the ratio value remained constant (k1/k2/E). Finally, the theoretical calculations were applied to the I26 backfill working face in the Xingdong mine to calculate the final subsidence amounts of the main roof. Field observations and theoretical calculations were about 48 mm, verifying the method’s applicability.


---

## Body

## 1. Introduction

Strata movement gradually develops on the surface after underground coal seam mining, causing surface subsidence, destruction of ground facilities, and so on. Solid backfill mining, paste backfill mining, high-water backfill mining, and other backfill mining methods are commonly used to prevent surface subsidence caused by underground mining. Qian et al. [1] put forward the concept of green mining in coal mine and expounded the basic methods of mining settlement and backfilling control. Wu et al. [2], through image algorithms and so on, analyzed the mining landscape changes before and after subsidence with the support of GIS technology. According to the characteristics of mining subsidence, Jung et al. studied the comprehensive prediction and calculation methods related to mining subsidence, and so on [3–9]. In view of the subsidence characteristics of the surface steps caused by the mining of the shallow and extra thick coal seam, Ju and Xu [10] put forward three possible control methods for surface stepped subsidence. Based on the field investigation and study of overburden damage in ultra-thick coal seam mining, the statistical formula was presented to estimate the maximum heights of failure zone in the LTCC operation [11]. The characteristics of displacement and ground subsidence caused by underground mining are studied by means of remote sensing and Geographic Information System [12–14]. Li et al. studied the mechanical characteristics of the gangue filling body, the mechanical and geometric characteristics of the hydraulic support for gangue filling, and the characteristics of the mining pressure in the solid filling face; they obtained that the solid backfilling method was an effective method to prevent hard-roof-induced face bursts, and the equivalent mining height model is capable of predicting surface deformation [15–20]. Benzaazoua et al. studied, respectively, the hardening process of cemented backfill, the preparation, microstructure characteristics, mechanical parameters, and backfilling technology of the paste backfill materials [21–24]. Based on the basic mechanical experiments of ultra-high-water materials, Ding et al. [25] found that it was good backfilling material. Yan et al. [26] put forward the pump filling technology of ZKD high-water quick-setting material, which was effectively applied to engineering practice. Feng et al. [27] introduced four kinds of backfilling methods for goaf with high water content materials and analyzed the advantages and disadvantages of each method.There have been many studies on the backfill mining process and issues that may arise. Zhang et al. [28] analyzed the interaction between the backfilling body and overburden strata in a fully mechanized backfilling mining face. Miao et al. [29] proposed the beam model on an elastic foundation with a fixed boundary condition to study the relationship between subsidence of the main roof and Young’s modulus of backfill body. Chen et al. [30] built a beam model on an elastic foundation with an elastic foundation boundary to analyze the variation of the main roof’s maximum subsidence with Young’s modulus of backfill body and elastic foundation coefficient of coal seam. Li et al. [31] proposed a thin plate model on an elastic foundation of the main roof with a fixed boundary condition to study the variation of the main roof’s subsidence and maximum tensile stress with an elastic foundation coefficient of the backfill body. Huang et al. [32] analyzed overlying strata movement in backfill mining using a similar physical simulation, and although these results are beneficial for scientific backfill mining, there are specific scenarios that they are applicable. The beam model in backfill mining is only suitable for roof mechanics analysis of the middle area in the working face with large length-width ratio.The thin plate model with the main roof being elastic and fixed boundary condition for backfill mining simplifies the main roof in order for it to be a fixed boundary condition. However, the main roof is bound to sink in the coal seam support area when the coal seam is thick and soft, making the fixed boundary condition unsuitable. The Winkler elastic foundation model should therefore be considered because of the weak shear-bearing capacity and non-deformability of the coal seam [33]. The main roof can then be regarded as an elastic foundation, and the elastic foundation model is able to analyze its mechanics and displacement characteristics during backfill mining.Here, we established a thin plate model on elastic foundation of the main roof with an elastic foundation boundary to study the variations of the internal force field and displacement field with various elastic foundation coefficients of the backfill body, elastic foundation coefficients of the coal seam, thicknesses of main roof, Young’s modulus of the main roof, and Poisson’s ratio of the main roof, and the fracture position and fracture conditions of the main roof were then determined and compared. Finally, a determination method for the elastic foundation coefficient of the backfill body is proposed in this paper. These theoretical calculations were then verified using the I26 mining face from Xingdong coal mine. These results have practical value in further developments in backfill mining technology.


## 2. The Thin Plate Model on Elastic Foundation with the Elastic Foundation Boundary in Backfill Mining

### 2.1. The Actual Surrounding Rock Condition of the Main Roof in Backfill Mining

The main roof is clamped between the overlying strata and immediate roof after mining. Due to the large compressive deformation and the weak shear-bearing capacity of the coal seam, the coal seam is best represented by the Winkler elastic foundation model. Similarly, the backfill body is also best represented by the Winkler elastic foundation model. Figure1 shows the surrounding rock relationship during backfill mining.Figure 1
Depiction of the surrounding rock in backfill mining.


### 2.2. Mechanical Hypothesis Condition of the Main Roof in Backfill Mining

#### 2.2.1. Hypothesis of Elastic Foundation Boundary

For this study, the coal seam is analyzed using the Winkler elastic foundation model. The displacement of the main roof is mainly limited by the immediate roof and coal seam. In general, Young’s modulus of the coal seam is much lower than that of the main roof and immediate roof; therefore, the coal seam is the main factor that limits the displacement of the main roof. The elastic foundation coefficientk0 can be derived from equation (1) and is given by(1)k0=E1E0E1h0+E0h1≈E1h1=k1,where k0 is the composite elastic foundation coefficient of the immediate roof and coal seam; k1 is the elastic foundation coefficient of the coal seam; h0 is the thickness of the immediate roof; E0 is Young’s modulus of the immediate roof; h1 is the thickness of the coal seam; E1 is Young’s modulus of the coal seam.


#### 2.2.2. Hypothesis of Elastic Foundation of the Backfill Body

Solid backfill material, paste backfill material, and high-water backfill material approximately satisfy the Winkler elastic foundation model.


#### 2.2.3. Dimension Requirement of the Elastic Thin Plate of the Main Roof

In order to satisfy the requirements to be considered a “thin plate,” the dimensions of the roof must meet the requirement of(2)18∼15≥hL≥1100∼180,where h is the thickness of the plate; L is the short side length of the plate.Whenh and L of the plate satisfy equation (2), the plate can be seen as a thin plate.In general, the ratio of the thickness of the roof to the short side length satisfies equation (2); that is, the main roof satisfies dimension requirements of the elastic thin plate.


### 2.3. Mechanical Model of the Thin Plate on Elastic Foundation with Elastic Foundation Boundary in Backfill Mining

#### 2.3.1. Mechanical Model

The mechanical model of a thin plate on an elastic foundation with an elastic foundation boundary of a main roof in backfill mining is shown in Figure2.Figure 2
The thin plate model on an elastic foundation with an elastic foundation boundary. (a) Vertical view; (b) I-I section plane; (c) II-II section plane.
(a)
(b)
(c)The rectangle area ABCD is the backfill areaS2. The length of AB is 2a and the length of AD is 2b. A coordinate system is established by using the midpoint of the main roof on the backfill body as the origin. Area S1 of the main roof, which is outside of ABCD and inside of A1B1C1D1, is approximately the clamping area of the elastic layer. The displacement of the main roof in area S1 is mainly determined by the stiffness of the coal seam. The load carried by the main roof is q0.


#### 2.3.2. Deflection Equation

According to the theory of a plate [34–36], the partial differential equation of the deflection of the main roof in an elastic foundation area S1 can be defined as(3)∂4ω1∂x4+2∂4ω1∂x2∂y2+∂4ω1∂y4=1D−k1ω1,where ω1 is the deflection function of the main roof in the area S1; and k1 is the elastic foundation coefficient of the coal seam. The stiffness of a thin plate is D = Eh3/(12−12μ2), where E is Young’s modulus of the main roof; h is the thickness of the main roof; and μ is Poisson’s ratio of the main roof.The partial differential equation of deflection of the main roof in the backfill areaS2 can be defined as(4)∂4ω2∂x4+2∂4ω2∂x2∂y2+∂4ω2∂y4=1Dq0−k2ω2,where ω2 is the deflection function of the main roof in area S2; q0 is the load carried by the main roof; and k2 is the elastic foundation coefficient of backfill body.


#### 2.3.3. Boundary Conditions

(1) Continuous Boundary Condition. Edges AB, BC, CD, and AD are the main roof’s interface between the coal body and the backfill body. Therefore, deflection, rotation angle, bending moment, and shear force are the same on the interface, and a continuous condition can be achieved by combining all the relevant equations as equations (5) and (6). A continuous boundary condition for the main roof can, therefore, be expressed as(5)−a≤x≤a,y=b,∂∂y∇2ω1=∂∂y∇2ω2,∂2ω1∂y2+μ∂2ω1∂x2=∂2ω2∂y2+μ∂2ω2∂x2,∂ω1∂y=∂ω2∂y,ω1x,b=ω2x,b,−a≤x≤a,y=−b,∂∂y∇2ω1=∂∂y∇2ω2,∂2ω1∂y2+μ∂2ω1∂x2=∂2ω2∂y2+μ∂2ω2∂x2,∂ω1∂y=∂ω2∂y,ω1x,−b=ω2x,−b,(6)−b≤y≤b,x=a,∂∂x∇2ω1=∂∂x∇2ω2,∂2ω1∂x2+μ∂2ω1∂y2=∂2ω2∂x2+μ∂2ω2∂y2,∂ω1∂x=∂ω2∂x,ω1a,y=ω2a,y,−b≤y≤b,x=−a,∂∂x∇2ω1=∂∂x∇2ω2,∂2ω1∂x2+μ∂2ω1∂y2=∂2ω2∂x2+μ∂2ω2∂y2,∂ω1∂x=∂ω2∂x,ω1−a,y=ω2−a,y.(2) Boundary Condition. A1B1 length is 2x0. A1D1 length is 2y0 (refer to Figure 2). At an infinite distance from the backfill area S2, the mining effect is very small. Therefore, the deflection and rotation angle of the main roof are both zero and satisfy the fixed boundary conditions at an infinite distance from the backfill area S2. The fixed boundary condition can then be defined as (7)x=−x0→−∞,−y0≤y≤y0,ω1=0,∂ω1∂x=0,x=x0→+∞,−y0≤y≤y0,ω1=0,∂ω1∂x=0,y=−y0→−∞,−x0≤x≤x0,ω1=0,∂ω1∂y=0,y=y0→+∞,−x0≤x≤x0,ω1=0,∂ω1∂y=0.


## 2.1. The Actual Surrounding Rock Condition of the Main Roof in Backfill Mining

The main roof is clamped between the overlying strata and immediate roof after mining. Due to the large compressive deformation and the weak shear-bearing capacity of the coal seam, the coal seam is best represented by the Winkler elastic foundation model. Similarly, the backfill body is also best represented by the Winkler elastic foundation model. Figure1 shows the surrounding rock relationship during backfill mining.Figure 1
Depiction of the surrounding rock in backfill mining.


## 2.2. Mechanical Hypothesis Condition of the Main Roof in Backfill Mining

### 2.2.1. Hypothesis of Elastic Foundation Boundary

For this study, the coal seam is analyzed using the Winkler elastic foundation model. The displacement of the main roof is mainly limited by the immediate roof and coal seam. In general, Young’s modulus of the coal seam is much lower than that of the main roof and immediate roof; therefore, the coal seam is the main factor that limits the displacement of the main roof. The elastic foundation coefficientk0 can be derived from equation (1) and is given by(1)k0=E1E0E1h0+E0h1≈E1h1=k1,where k0 is the composite elastic foundation coefficient of the immediate roof and coal seam; k1 is the elastic foundation coefficient of the coal seam; h0 is the thickness of the immediate roof; E0 is Young’s modulus of the immediate roof; h1 is the thickness of the coal seam; E1 is Young’s modulus of the coal seam.


### 2.2.2. Hypothesis of Elastic Foundation of the Backfill Body

Solid backfill material, paste backfill material, and high-water backfill material approximately satisfy the Winkler elastic foundation model.


### 2.2.3. Dimension Requirement of the Elastic Thin Plate of the Main Roof

In order to satisfy the requirements to be considered a “thin plate,” the dimensions of the roof must meet the requirement of(2)18∼15≥hL≥1100∼180,where h is the thickness of the plate; L is the short side length of the plate.Whenh and L of the plate satisfy equation (2), the plate can be seen as a thin plate.In general, the ratio of the thickness of the roof to the short side length satisfies equation (2); that is, the main roof satisfies dimension requirements of the elastic thin plate.


## 2.2.1. Hypothesis of Elastic Foundation Boundary

For this study, the coal seam is analyzed using the Winkler elastic foundation model. The displacement of the main roof is mainly limited by the immediate roof and coal seam. In general, Young’s modulus of the coal seam is much lower than that of the main roof and immediate roof; therefore, the coal seam is the main factor that limits the displacement of the main roof. The elastic foundation coefficientk0 can be derived from equation (1) and is given by(1)k0=E1E0E1h0+E0h1≈E1h1=k1,where k0 is the composite elastic foundation coefficient of the immediate roof and coal seam; k1 is the elastic foundation coefficient of the coal seam; h0 is the thickness of the immediate roof; E0 is Young’s modulus of the immediate roof; h1 is the thickness of the coal seam; E1 is Young’s modulus of the coal seam.


## 2.2.2. Hypothesis of Elastic Foundation of the Backfill Body

Solid backfill material, paste backfill material, and high-water backfill material approximately satisfy the Winkler elastic foundation model.


## 2.2.3. Dimension Requirement of the Elastic Thin Plate of the Main Roof

In order to satisfy the requirements to be considered a “thin plate,” the dimensions of the roof must meet the requirement of(2)18∼15≥hL≥1100∼180,where h is the thickness of the plate; L is the short side length of the plate.Whenh and L of the plate satisfy equation (2), the plate can be seen as a thin plate.In general, the ratio of the thickness of the roof to the short side length satisfies equation (2); that is, the main roof satisfies dimension requirements of the elastic thin plate.


## 2.3. Mechanical Model of the Thin Plate on Elastic Foundation with Elastic Foundation Boundary in Backfill Mining

### 2.3.1. Mechanical Model

The mechanical model of a thin plate on an elastic foundation with an elastic foundation boundary of a main roof in backfill mining is shown in Figure2.Figure 2
The thin plate model on an elastic foundation with an elastic foundation boundary. (a) Vertical view; (b) I-I section plane; (c) II-II section plane.
(a)
(b)
(c)The rectangle area ABCD is the backfill areaS2. The length of AB is 2a and the length of AD is 2b. A coordinate system is established by using the midpoint of the main roof on the backfill body as the origin. Area S1 of the main roof, which is outside of ABCD and inside of A1B1C1D1, is approximately the clamping area of the elastic layer. The displacement of the main roof in area S1 is mainly determined by the stiffness of the coal seam. The load carried by the main roof is q0.


### 2.3.2. Deflection Equation

According to the theory of a plate [34–36], the partial differential equation of the deflection of the main roof in an elastic foundation area S1 can be defined as(3)∂4ω1∂x4+2∂4ω1∂x2∂y2+∂4ω1∂y4=1D−k1ω1,where ω1 is the deflection function of the main roof in the area S1; and k1 is the elastic foundation coefficient of the coal seam. The stiffness of a thin plate is D = Eh3/(12−12μ2), where E is Young’s modulus of the main roof; h is the thickness of the main roof; and μ is Poisson’s ratio of the main roof.The partial differential equation of deflection of the main roof in the backfill areaS2 can be defined as(4)∂4ω2∂x4+2∂4ω2∂x2∂y2+∂4ω2∂y4=1Dq0−k2ω2,where ω2 is the deflection function of the main roof in area S2; q0 is the load carried by the main roof; and k2 is the elastic foundation coefficient of backfill body.


### 2.3.3. Boundary Conditions

(1) Continuous Boundary Condition. Edges AB, BC, CD, and AD are the main roof’s interface between the coal body and the backfill body. Therefore, deflection, rotation angle, bending moment, and shear force are the same on the interface, and a continuous condition can be achieved by combining all the relevant equations as equations (5) and (6). A continuous boundary condition for the main roof can, therefore, be expressed as(5)−a≤x≤a,y=b,∂∂y∇2ω1=∂∂y∇2ω2,∂2ω1∂y2+μ∂2ω1∂x2=∂2ω2∂y2+μ∂2ω2∂x2,∂ω1∂y=∂ω2∂y,ω1x,b=ω2x,b,−a≤x≤a,y=−b,∂∂y∇2ω1=∂∂y∇2ω2,∂2ω1∂y2+μ∂2ω1∂x2=∂2ω2∂y2+μ∂2ω2∂x2,∂ω1∂y=∂ω2∂y,ω1x,−b=ω2x,−b,(6)−b≤y≤b,x=a,∂∂x∇2ω1=∂∂x∇2ω2,∂2ω1∂x2+μ∂2ω1∂y2=∂2ω2∂x2+μ∂2ω2∂y2,∂ω1∂x=∂ω2∂x,ω1a,y=ω2a,y,−b≤y≤b,x=−a,∂∂x∇2ω1=∂∂x∇2ω2,∂2ω1∂x2+μ∂2ω1∂y2=∂2ω2∂x2+μ∂2ω2∂y2,∂ω1∂x=∂ω2∂x,ω1−a,y=ω2−a,y.(2) Boundary Condition. A1B1 length is 2x0. A1D1 length is 2y0 (refer to Figure 2). At an infinite distance from the backfill area S2, the mining effect is very small. Therefore, the deflection and rotation angle of the main roof are both zero and satisfy the fixed boundary conditions at an infinite distance from the backfill area S2. The fixed boundary condition can then be defined as (7)x=−x0→−∞,−y0≤y≤y0,ω1=0,∂ω1∂x=0,x=x0→+∞,−y0≤y≤y0,ω1=0,∂ω1∂x=0,y=−y0→−∞,−x0≤x≤x0,ω1=0,∂ω1∂y=0,y=y0→+∞,−x0≤x≤x0,ω1=0,∂ω1∂y=0.


## 2.3.1. Mechanical Model

The mechanical model of a thin plate on an elastic foundation with an elastic foundation boundary of a main roof in backfill mining is shown in Figure2.Figure 2
The thin plate model on an elastic foundation with an elastic foundation boundary. (a) Vertical view; (b) I-I section plane; (c) II-II section plane.
(a)
(b)
(c)The rectangle area ABCD is the backfill areaS2. The length of AB is 2a and the length of AD is 2b. A coordinate system is established by using the midpoint of the main roof on the backfill body as the origin. Area S1 of the main roof, which is outside of ABCD and inside of A1B1C1D1, is approximately the clamping area of the elastic layer. The displacement of the main roof in area S1 is mainly determined by the stiffness of the coal seam. The load carried by the main roof is q0.


## 2.3.2. Deflection Equation

According to the theory of a plate [34–36], the partial differential equation of the deflection of the main roof in an elastic foundation area S1 can be defined as(3)∂4ω1∂x4+2∂4ω1∂x2∂y2+∂4ω1∂y4=1D−k1ω1,where ω1 is the deflection function of the main roof in the area S1; and k1 is the elastic foundation coefficient of the coal seam. The stiffness of a thin plate is D = Eh3/(12−12μ2), where E is Young’s modulus of the main roof; h is the thickness of the main roof; and μ is Poisson’s ratio of the main roof.The partial differential equation of deflection of the main roof in the backfill areaS2 can be defined as(4)∂4ω2∂x4+2∂4ω2∂x2∂y2+∂4ω2∂y4=1Dq0−k2ω2,where ω2 is the deflection function of the main roof in area S2; q0 is the load carried by the main roof; and k2 is the elastic foundation coefficient of backfill body.


## 2.3.3. Boundary Conditions

(1) Continuous Boundary Condition. Edges AB, BC, CD, and AD are the main roof’s interface between the coal body and the backfill body. Therefore, deflection, rotation angle, bending moment, and shear force are the same on the interface, and a continuous condition can be achieved by combining all the relevant equations as equations (5) and (6). A continuous boundary condition for the main roof can, therefore, be expressed as(5)−a≤x≤a,y=b,∂∂y∇2ω1=∂∂y∇2ω2,∂2ω1∂y2+μ∂2ω1∂x2=∂2ω2∂y2+μ∂2ω2∂x2,∂ω1∂y=∂ω2∂y,ω1x,b=ω2x,b,−a≤x≤a,y=−b,∂∂y∇2ω1=∂∂y∇2ω2,∂2ω1∂y2+μ∂2ω1∂x2=∂2ω2∂y2+μ∂2ω2∂x2,∂ω1∂y=∂ω2∂y,ω1x,−b=ω2x,−b,(6)−b≤y≤b,x=a,∂∂x∇2ω1=∂∂x∇2ω2,∂2ω1∂x2+μ∂2ω1∂y2=∂2ω2∂x2+μ∂2ω2∂y2,∂ω1∂x=∂ω2∂x,ω1a,y=ω2a,y,−b≤y≤b,x=−a,∂∂x∇2ω1=∂∂x∇2ω2,∂2ω1∂x2+μ∂2ω1∂y2=∂2ω2∂x2+μ∂2ω2∂y2,∂ω1∂x=∂ω2∂x,ω1−a,y=ω2−a,y.(2) Boundary Condition. A1B1 length is 2x0. A1D1 length is 2y0 (refer to Figure 2). At an infinite distance from the backfill area S2, the mining effect is very small. Therefore, the deflection and rotation angle of the main roof are both zero and satisfy the fixed boundary conditions at an infinite distance from the backfill area S2. The fixed boundary condition can then be defined as (7)x=−x0→−∞,−y0≤y≤y0,ω1=0,∂ω1∂x=0,x=x0→+∞,−y0≤y≤y0,ω1=0,∂ω1∂x=0,y=−y0→−∞,−x0≤x≤x0,ω1=0,∂ω1∂y=0,y=y0→+∞,−x0≤x≤x0,ω1=0,∂ω1∂y=0.


## 3. The Finite Difference Method for Solving the Partial Differential Equations

As mentioned above, it is extremely difficult to obtain the exact solutions of the above partial differential equations, but the approximate solutions satisfy engineering and practical mining requirement. The finite difference method is an effective method to obtain the approximate solutions of differential equations [37–39].


### 3.1. The Finite Difference Method

#### 3.1.1. The Nodal Layout of the Difference Equation

According to the finite difference theory, the difference equation involving 13 nodes is needed to solve the partial differential equation [34].The difference grid layout of 13 nodes is shown in Figure3. Δx = Δy = d is the nodal space. Point J0 is a feature node. Deflection of point J0 is expressed as ωij. Numbers of the remaining nodes are determined by the intersection of the vertical and horizontal lines (refer to Figure 3).Figure 3
The difference grid layout of 13 nodes.


#### 3.1.2. The Difference Equations of the Partial Differential Equations

Combining with the difference grid layout and the finite difference method of the partial differential equations (refer to Figure3), the difference equations of the partial differential equations (3) and (4) at feature node J0 can be obtained by(8)20+d4k1Dωij−8ωi+1,j+ωi−1,j+ωi,j+1+ωi,j−1+2ωi+1,j+1+ωi+1,j−1+ωi−1,j+1+ωi−1,j−1+ωi+2,j+ωi−2,j+ωi,j+2+ωi,j−2=0,where d is the nodal space; ωij, ωi+1,j, ωi−1,j, ωi,j+1, ωi,j−1, ωi+1,j+1, ωi+1,j−1, ωi−1,j+1, ωi−1,j−1, ωi+2,j, ωi−2,j, ωi,j+2 and ωi,j−2 are nodal deflection, respectively (refer to Figure 3).(9)20+d4k2Dωij−8ωi+1+ωi−1,j+ωi,j+1+ωi,j−1+2ωi+1j+1+ωi+1,j−1+ωi−1,j+1+ωi−1,j−1+ωi+2,j+ωi−2,j+ωi,j+2+ωi,j−2=q0d4D.


#### 3.1.3. The Difference Equations of Boundary Condition

Theoretically, the exact fixed boundary condition can be satisfied whenx0 and y0 approach infinity (refer to equation (7)), but the finite difference method cannot calculate this infinite region. Actually, when the outer boundary is three times the maximum length of the mining area (y0 = x0 = 3max{2a, 2b}), the influence of mining on this boundary is very weak, which approximately satisfies both the fixed boundary condition and engineering requirements. Nodal space d is taken as 0.2 m. The difference equations of the boundary condition in feature node J0 can be expressed as(10)−y0≤y≤y0,x=−x0,ωi,j=0,∂ω∂xij=ωi−1,j−ωi+1,j2d=0,−y0≤y≤y0,x=x0,ωi,j=0,∂ω∂xij=ωi−1,j−ωi+1,j2d=0,−x0≤x≤x0,y=y0,ωi,j=0,∂ω∂yij=ωi,j−1−ωi,j+12d=0,−x0≤x≤x0,y=−y0,ωi,j=0,∂ω∂yij=ωi,j−1−ωi,j+12d=0.


### 3.2. Calculation Process

According to the finite difference method, the 13-node difference equation should be established for any node whose deflection is unknown within the area ofS1 and the S2 area of the main roof. The nodes whose deflection is unknown can then be obtained by constructing all nodal difference equations and the boundary condition equations.Figure4 presents the specific process of solving the partial differential equations. The process of solving partial differential equations needs to adopt the Sparse function in Matlab to form algebraic equations for a sparse matrix, and each nodal deflection solution can be obtained by using function Gmres [37] to solve the algebraic equations.Figure 4
Process of solving partial differential equations.After solving each nodal deflection of the main roof, the nodal internal force solutions can be obtained by(11)Mxij=−D∂2ω∂x2+μ∂2ω∂y2ij=−Dd2ωi−1,j−2ωi,j+ωi+1,j−μωi,j−1−2ωi,j+ωi,j+1,Myij=−D∂2ω∂y2+μ∂2ω∂x2ij=−Dd2ωi,j−1−2ωi,j+ωi,j+1−μωi−1,j−2ωi,j+ωi+1,j,Mxyij=−D1−μ∂2ω∂x∂yij=−D1−μ4d2ωi−1,j−1−ωi+1,−1+ωi+1,j+1−ωi−1,j+1,where (Mx)ij is the nodal moment of the x-component; (My)ij is the nodal moment of the y-component; and (Mxy)ij is the nodal twisting moment.


## 3.1. The Finite Difference Method

### 3.1.1. The Nodal Layout of the Difference Equation

According to the finite difference theory, the difference equation involving 13 nodes is needed to solve the partial differential equation [34].The difference grid layout of 13 nodes is shown in Figure3. Δx = Δy = d is the nodal space. Point J0 is a feature node. Deflection of point J0 is expressed as ωij. Numbers of the remaining nodes are determined by the intersection of the vertical and horizontal lines (refer to Figure 3).Figure 3
The difference grid layout of 13 nodes.


### 3.1.2. The Difference Equations of the Partial Differential Equations

Combining with the difference grid layout and the finite difference method of the partial differential equations (refer to Figure3), the difference equations of the partial differential equations (3) and (4) at feature node J0 can be obtained by(8)20+d4k1Dωij−8ωi+1,j+ωi−1,j+ωi,j+1+ωi,j−1+2ωi+1,j+1+ωi+1,j−1+ωi−1,j+1+ωi−1,j−1+ωi+2,j+ωi−2,j+ωi,j+2+ωi,j−2=0,where d is the nodal space; ωij, ωi+1,j, ωi−1,j, ωi,j+1, ωi,j−1, ωi+1,j+1, ωi+1,j−1, ωi−1,j+1, ωi−1,j−1, ωi+2,j, ωi−2,j, ωi,j+2 and ωi,j−2 are nodal deflection, respectively (refer to Figure 3).(9)20+d4k2Dωij−8ωi+1+ωi−1,j+ωi,j+1+ωi,j−1+2ωi+1j+1+ωi+1,j−1+ωi−1,j+1+ωi−1,j−1+ωi+2,j+ωi−2,j+ωi,j+2+ωi,j−2=q0d4D.


### 3.1.3. The Difference Equations of Boundary Condition

Theoretically, the exact fixed boundary condition can be satisfied whenx0 and y0 approach infinity (refer to equation (7)), but the finite difference method cannot calculate this infinite region. Actually, when the outer boundary is three times the maximum length of the mining area (y0 = x0 = 3max{2a, 2b}), the influence of mining on this boundary is very weak, which approximately satisfies both the fixed boundary condition and engineering requirements. Nodal space d is taken as 0.2 m. The difference equations of the boundary condition in feature node J0 can be expressed as(10)−y0≤y≤y0,x=−x0,ωi,j=0,∂ω∂xij=ωi−1,j−ωi+1,j2d=0,−y0≤y≤y0,x=x0,ωi,j=0,∂ω∂xij=ωi−1,j−ωi+1,j2d=0,−x0≤x≤x0,y=y0,ωi,j=0,∂ω∂yij=ωi,j−1−ωi,j+12d=0,−x0≤x≤x0,y=−y0,ωi,j=0,∂ω∂yij=ωi,j−1−ωi,j+12d=0.


## 3.1.1. The Nodal Layout of the Difference Equation

According to the finite difference theory, the difference equation involving 13 nodes is needed to solve the partial differential equation [34].The difference grid layout of 13 nodes is shown in Figure3. Δx = Δy = d is the nodal space. Point J0 is a feature node. Deflection of point J0 is expressed as ωij. Numbers of the remaining nodes are determined by the intersection of the vertical and horizontal lines (refer to Figure 3).Figure 3
The difference grid layout of 13 nodes.


## 3.1.2. The Difference Equations of the Partial Differential Equations

Combining with the difference grid layout and the finite difference method of the partial differential equations (refer to Figure3), the difference equations of the partial differential equations (3) and (4) at feature node J0 can be obtained by(8)20+d4k1Dωij−8ωi+1,j+ωi−1,j+ωi,j+1+ωi,j−1+2ωi+1,j+1+ωi+1,j−1+ωi−1,j+1+ωi−1,j−1+ωi+2,j+ωi−2,j+ωi,j+2+ωi,j−2=0,where d is the nodal space; ωij, ωi+1,j, ωi−1,j, ωi,j+1, ωi,j−1, ωi+1,j+1, ωi+1,j−1, ωi−1,j+1, ωi−1,j−1, ωi+2,j, ωi−2,j, ωi,j+2 and ωi,j−2 are nodal deflection, respectively (refer to Figure 3).(9)20+d4k2Dωij−8ωi+1+ωi−1,j+ωi,j+1+ωi,j−1+2ωi+1j+1+ωi+1,j−1+ωi−1,j+1+ωi−1,j−1+ωi+2,j+ωi−2,j+ωi,j+2+ωi,j−2=q0d4D.


## 3.1.3. The Difference Equations of Boundary Condition

Theoretically, the exact fixed boundary condition can be satisfied whenx0 and y0 approach infinity (refer to equation (7)), but the finite difference method cannot calculate this infinite region. Actually, when the outer boundary is three times the maximum length of the mining area (y0 = x0 = 3max{2a, 2b}), the influence of mining on this boundary is very weak, which approximately satisfies both the fixed boundary condition and engineering requirements. Nodal space d is taken as 0.2 m. The difference equations of the boundary condition in feature node J0 can be expressed as(10)−y0≤y≤y0,x=−x0,ωi,j=0,∂ω∂xij=ωi−1,j−ωi+1,j2d=0,−y0≤y≤y0,x=x0,ωi,j=0,∂ω∂xij=ωi−1,j−ωi+1,j2d=0,−x0≤x≤x0,y=y0,ωi,j=0,∂ω∂yij=ωi,j−1−ωi,j+12d=0,−x0≤x≤x0,y=−y0,ωi,j=0,∂ω∂yij=ωi,j−1−ωi,j+12d=0.


## 3.2. Calculation Process

According to the finite difference method, the 13-node difference equation should be established for any node whose deflection is unknown within the area ofS1 and the S2 area of the main roof. The nodes whose deflection is unknown can then be obtained by constructing all nodal difference equations and the boundary condition equations.Figure4 presents the specific process of solving the partial differential equations. The process of solving partial differential equations needs to adopt the Sparse function in Matlab to form algebraic equations for a sparse matrix, and each nodal deflection solution can be obtained by using function Gmres [37] to solve the algebraic equations.Figure 4
Process of solving partial differential equations.After solving each nodal deflection of the main roof, the nodal internal force solutions can be obtained by(11)Mxij=−D∂2ω∂x2+μ∂2ω∂y2ij=−Dd2ωi−1,j−2ωi,j+ωi+1,j−μωi,j−1−2ωi,j+ωi,j+1,Myij=−D∂2ω∂y2+μ∂2ω∂x2ij=−Dd2ωi,j−1−2ωi,j+ωi,j+1−μωi−1,j−2ωi,j+ωi+1,j,Mxyij=−D1−μ∂2ω∂x∂yij=−D1−μ4d2ωi−1,j−1−ωi+1,−1+ωi+1,j+1−ωi−1,j+1,where (Mx)ij is the nodal moment of the x-component; (My)ij is the nodal moment of the y-component; and (Mxy)ij is the nodal twisting moment.


## 4. Example and Analysis

The thin plate model on an elastic foundation with an elastic foundation boundary in backfill mining was investigated under realistic geological conditions. Variations of the internal force and deflection of the main roof withE, h, μ, k1 and k2 were studied using the above equations for theoretical calculations.The I26 working face of the Xingdong mine (length AB = 70 m and advancing length AD = 200 m in Figure2) adopted the backfill bag of high-water material to backfill mining, with a backfill rate of 100%. This working face was used to obtain realistic conditions for the calculations. Based on the above analysis, when A1B1 = A1D1 = 600 m, the outer boundary can be regarded as fixed boundary.In the I26 working face, the average angle of the coal seam was 4–6°, mining height was 4.5 m, and the thickness of the main roof was 10 m. Young’s modulus of the main roof was 26 Gpa and the elastic foundation coefficientk1 was 1000 MN/m3. The key rock thickness, which overlies the main roof, was 110 m and the average unit weight of the rock was 0.024 MN/m3. Therefore, the load carried by the main roof was q0 = 2.64 Mpa. The elastic foundation coefficient k2 of the high-water material used during the backfill mining process was 20 MN/m3–220 MN/m3.


### 4.1. Fracture Criterion

(12)M1,M3=Mx+My2±Mx−My22+Mxy2,where M1 is the maximum principal bending moment; M3 is the minimum principal bending moment; Mx is the moment of x-component; My is the moment of y-component; and Mxy is the nodal twisting moment.(13)M1ijM3ij=Mxij+Myij2±Mxij−Myij22+Mxyij2,where (M1)ij is the nodal maximum principal moment; (M3)ij is the nodal minimum principal moment.Because the tensile strength of the rock is much less than the compressive strength, the principal moment of rock is compared with the ultimate bending moment of rock to judge whether the rock has failed or not. The maximum principal momentM1 and minimum principal moment M3 are obtained by equation (12), and the difference equation of equation (12) is equation (13). The moments (Mx)ij, (My)ij and (Mxy)ij in equation (13) can be obtained by equation (11). The extreme points of the principal moment are extracted from the cloud chart of the principal moment, which are shown in Figure 5.Figure 5
Principal moment feature whenk2 = 20 MN/m3. (a) Three-dimensional cloud chart of principal moment M1; (b) plane cloud chart of principal moment M1; (c) cloud chart of principal moment M3; (d) plane cloud chart of principal moment M3.
(a)
(b)
(c)
(d)


### 4.2. Effect of Elastic Foundation Coefficientk2 of Backfill Body in Backfill Mining

#### 4.2.1. Basic Principles and Distribution Pattern of the Principal Moment

The elastic foundation coefficient of the backfill body is the key material factor to decide the backfill effect. Here, the backfill rate was 100%. The elastic foundation coefficientk2 was ∼20 MN/m3–220 MN/m3. When k2 = 20 MN/m3 and k2 = 180 MN/m3, cloud charts of the principal moment in backfill mining were plotted, as shown in Figures 5 and 6.Figure 6
Principal moment feature whenk2 = 180 MN/m3. (a) Three-dimensional cloud chart of principal moment M1; (b) plane cloud chart of principal moment M1; (c) three-dimensional cloud chart of principal moment M3; (d) plane cloud chart of principal moment M3.
(a)
(b)
(c)
(d)From Figures5 and 6, we can obtain the distribution of the principal moment during backfill mining: principal moments M1 and M3 in the middle area of the main roof directly overlying the backfill body were all positive bending moments; namely, the lower surface of the main roof bores the tensile stress and the upper surface of the main roof bores the compressive stress in this area. The nodal principal moments M3 on the long side and short side area with a certain distance to the coal wall (Lc) were negative bending moments; namely, the upper surface of the main roof bores the tensile stress and the lower surface of main roof bores the compressive stress in this area.According to the distribution characteristics of the principal moment in Figures5 and 6, the coordinates of the extreme points of the main roof’s principal moment were plotted, as shown in Figure 7.Figure 7
The coordinate of the extreme points of the main roof principal moment during backfill mining: (a)k2 = 20 MN/m3; (b) k2 = 180 MN/m3.
(a)
(b)When the value ofk2 was smaller (k2 = 20 MN/m3), and according to the distribution characteristics of the principal moment in Figure 5, the coordinates of the extreme points of the principal moment were plotted, as shown in Figure 7(a). The maximum absolute value of the principal moment of main roof in backfill area was at the midpoint with a coordinate of (0, 0), which was set to Mz, and Mz = M1|(0, 0); the maximum absolute value of the principal moment on the long side area with a certain distance to the coal wall (Lc) was the negative value of the minimum principal moment M3, which was set to Mc, the coordinates of which were (0, b + Lc) and (0, −b−Lc). The distance of the Mc advancing coal wall was set to Lc, and Mc=−M30,b+Lc=−M30,−b−Lc; the maximum absolute value of the principal moment on the short side area with a certain distance to the coal wall (Lc) was the negative value of the minimum principal moment M3, which was set to Md, and the coordinates of Md are (a + Ld, 0) and (−a−Ld, 0). The distance of Md advancing coal wall was set to Ld, and Md=−M3a+Ld,0=−M3−a−Ld,0.The coordinates of the extreme points of the principal moment from Figure6 were plotted when k2 was larger (k2 = 180 MN/m3), as shown in Figure 7(b). The maximum principal moment of main roof in backfill area was set to Mzc whose coordinates were (−a + Ldc, b−Lcc) and because of the symmetry, there were four points whose moments were all taken as MzcMzc−a+Ldc,b−Lcc>Mz=M10,0. In Figure 7(b), the maximum absolute value of the principal moment of the long side area with a certain distance to the coal wall (Lc) was set to Mc, and Mc=−M30,b+Lc=−M30,−b−Lc; the maximum absolute value of the principal moment of the short side area with a certain distance to the coal wall (Lc) wall was set to Md (Md=−M3a+Ld,0=−M3−b−Lc,0).As shown in Figures5 and 6, when the value of k2 became large (k2 = 180 MN/m3), the maximum principal moment Mz in the midpoint of the main roof transferred to the surrounding area and the points of the maximum principal moment were in the four corners of the study area. The maximum absolute value of the principal moment of main roof in the backfill area is denoted as Mzz, and Mzz meets the following relationship: Mzz=maxM10,0,Mzc−a+Ldc,b−LccBased on the above analysis, there appear to be three types of extreme points of the principal moments: the maximum absolute value of the principal momentMzz of the main roof in backfill area, the maximum absolute value of the principal moment Mc in the long side area, and the maximum absolute value of the principal moment Md in the short side area. Based on the above analysis and rock tensile properties, the lower surface of the middle area of main roof and the upper surface of long side and short side areas in front of the coal wall are the most susceptible to fracturing. Using the plate model on an elastic foundation with an elastic foundation boundary during backfill mining, the principles governing the fracturing of main roof are determined by the values of Mc, Mzz and Md. The main roof can only fracture when the principal moment is larger than Ms (the ultimate bending moment of the main roof).During backfill mining, the principal moments and the maximum subsidenceω0 of the main roof are related to k2, k1, h, E and μ. The principal movements and Lc were then calculated with variations in the above mentioned variables and analyzed in order to study the variations in the maximum subsidence of the main roof.


#### 4.2.2. Effects of the Elastic Foundation Coefficientk2 of the Backfill Body

The following calculation examines the effect of the elastic foundation coefficientk2 on the principal bending moment and maximum subsidence ω0 of the main roof.(1)
The variation of the principal bending moment of the main roof andLc with the elastic foundation coefficient k2 of the backfill bodyFigure8 shows the variations of the principal bending moments of the main roof and Lc with the elastic foundation coefficient k2 of the backfill body. The backfill body significantly reduced the value of the maximum principal bending moment of the main roof and improved the stability of the main roof. The principal bending moments (Mc, Md and Mzz) continuously decreased with the increase in values of k2. The value of the principal bending moment Mzz was highly sensitive to k2, and comparatively the principal bending moment Md was less sensitive. The principal bending moments Mc and Md were finally almost equal with an increase in k2. In this situation, the long side and short side of the main roof almost fractured simultaneously when the maximum principal bending moment of the main roof increased to the ultimate bending moment Ms. There should be very little decrease in the principal bending moments when k2 increases to some extent.The initial fracture position of the main roof was in the long side area, which is located in front of the coal wall with aMc value greater than Ms. Lc continuously increased with the increase in k2. Specifically, when Mc was greater than Ms, the larger the values of k2, the further the distance between the fracture line of main roof and the coal wall (refer to Figure 8).(2)
Variation of the maximum subsidenceω0 of the main roof with an elastic foundation coefficient k2 of the backfill bodyFigure9 shows the variation of the maximum subsidence ω0 of the main roof with an elastic foundation coefficient k2 of the backfill body. The maximum subsidence ω0 of the main roof appeared to decrease with the increase in k2. The maximum subsidence should decrease very little when k2 increases to some extent.Using the above calculations and formulas, further study should be conducted to further describe the principal bending moments and the maximum subsidenceω0 of main roof with an elastic foundation coefficient k1 by varying the coal seam, thickness h of the main roof, Young’s modulus E of the main roof, and Poisson’s ratio μ of the main roof, in order to better and widely study fracture principles and subsidence of the main roof during backfill mining.Figure 8
Change in the principal bending moments andLc with elastic foundation coefficient k2 of backfill body.Figure 9
Change in the maximum subsidenceω0 of the main roof with variable k2.


### 4.3. Effect of the Elastic Foundation Coefficientk1 of the Coal Seam in Backfill Mining

The elastic foundation coefficientk1 of the coal seam mainly reflects the limitation of the displacement of the main roof. We then analyzed variations of the principal bending moments and the maximum subsidence ω0 with an elastic foundation coefficient k1.(1)
Variations in the principal bending moments andLc with k1.Variations in the principal bending moments of the main roof and the distance to the coal wallLc with an elastic foundation coefficient k1 are presented in Figure 10. The principal bending moments (Mc, Md and Mzz) apparently increased with the increase in k1, and comparatively the principal bending moment Mzz had a smaller growth rate. When the principal moment Mc was greater than the ultimate moment Ms, the larger the value of k1 (the harder the coal seam), the closer the distance between the fracture line of the main roof and the coal wall (refer to Figure 10). Therefore, the influence of the boundary condition on the principal bending moment of main roof and Lc cannot be neglected. In other words, the influence of the boundary condition on the fracture position of main roof cannot be neglected.(2)
Variations in the maximum subsidenceω0 of the main roof with an elastic foundation coefficient k1 of coal seam.The change in the maximum subsidenceω0 with an elastic foundation coefficient k1 of the coal seam is presented in Figure 11. The maximum subsidence continuously decreased with an increase in k1, whereas the sensitivity of the maximum subsidence ω0 decreased with an increase in k1. The value of the maximum subsidence was reduced by approximately 10% with an increase in k1 (varying from 1000 to 8000 MN/m3).Figure 10
Change in the principal bending moments andLc with k1.Figure 11
Change in the maximum subsidenceω0 with k1.


### 4.4. Effect of the Thicknessh of the Main Roof

The thickness of the main roof can vary from mine to mine, depending on various circumstances (geology, previous mining, etc.). The effects of various thicknessesh were analyzed on the maximum subsidence ω0 and the principal bending moments in order to encompass a wide range of scenarios.(1)
Change in the principal bending moments andLc with hFigure12 graphically presents the variations of each of the three principal bending moments and Lc with various thicknesses of the main roof h (4–12 m). The principal bending moments (Mc, Md and Mzz) continuously increased with the increase in the roof thickness. Principal bending moment Mzz was highly sensitive to h, and comparatively the principal bending moment Md was less sensitive. The principal bending moments Mc and Md were almost equal when the thickness was relatively small (4–8 m). In such situation, the long side and short side area of the main roof almost fractured concurrently when the maximum principal bending moment of main roof increased to the ultimate bending moment Ms. Lc increased linearly with the increase in the roof thickness. Actually, the thicker the roof, the further the distance between the fracture line of main roof and the coal wall (refer to Figure 12).(2)
Changes in the maximum subsidence of the main roof with various thicknesseshFigure13 shows the change in the maximum subsidence ω0 of the main roof with various thicknesses (4–12 m). The maximum subsidence of the main roof continuously decreased with an increase in roof thickness, and the sensitivity of the maximum subsidence increased with an increase in the roof thickness.Figure 12
Change in the principal bending moments andLc with h.Figure 13
The maximum subsidenceω0 of main roof with roof thickness h.


### 4.5. Variations of Young’s ModulusE of the Main Roof in Backfill Mining

Young’s modulusE of the main roof was then varied to analyze the effects on the principal bending moment and maximum subsidence ω0.(1)
Changes in the principal bending moments andLc with various EFigure14 shows the changes in the principal bending moments and Lc with various Young’s modulus values E (10–75 GPa). The principal bending moments Mc and Md were almost equal when E was relatively small (10 GPa). In this condition, the long side and short side area of the main roof almost fractured concurrently when the maximum principal bending moment of the main roof increased to Ms. When the value of Young’s modulus reached 70 GPa, Mc and Mzz were nearly equal. In such situation, the long side and the middle area of the main roof almost fractured simultaneously when the principal bending moment was greater than Ms. The larger the E, the farther the Lc, and the farther the distance between the fracture line of main roof and the coal wall (refer to Figure 14). Based on the effect of Young’s modulus on the main bending moment, we can conclude that the greater Young’s modulus E value and the tensile strength must also be large within the main roof to ensure the main roof does not fracture.(2)
Change in maximum subsidenceω0 with various Young’s modulus E values of the main roofFigure15 shows the change in the maximum subsidence ω0 with an increase in Young’s modulus E of the main roof. Maximum subsidence of the main roof continuously decreased with an increase in E, and the sensitivity of ω0 decreased with an increase in E.Figure 14
Change in the principal bending moments andLc with E.Figure 15
The change in the maximum subsidenceω0 with various E.


### 4.6. Effects of Poisson’s Ratioμ on the Main Roof in Backfill Mining

Variations of the principal moments and maximum subsidenceω0 with Poisson’s ratio μ of main roof are discussed in the following calculations. Generally, the value of Poisson’s ratio varies from 0.1 to 0.3.(1)
Changes in the principal moments andLc with various Poisson’s ratio μ valuesFigure16 shows variations in principal moments and Lc with changes in Poisson’s ratio μ (0.10–0.30). The principal moments (Mc, Md and Mzz) gradually increased with the increase in μ, and Lc remained almost unchanged. Poisson’s ratio μ had little effect on principal moments and Lc; therefore, Poisson’s ratio would have little effect on the position of fracturing in the main roof.(2)
Variations of maximum subsidenceω0 of the main roof with μFigure17 shows how the maximum subsidence ω0 of the main roof varies with μ (0.10–0.30). The maximum subsidence slowly decreased with an increase in Poisson’s ratio. The maximum subsidence values ω0 when μ = 0.1 and μ = 0.3 were almost equal. Overall, Poisson’s ratio has little effect on the maximum subsidence of the main roof.Figure 16
The variation of principal bending moments andLc with μ.Figure 17
Variations of the maximum subsidenceω0 with μ.


### 4.7. Effects ofk1, k2, E in Backfill Mining

As observed through the above calculations, both the distribution of the principal bending moments of the main roof and the distance of the advancing coal seamLc were influenced by k1, k2 and E. The value of the principal bending moments decreased as k2 increased while the principal bending moments increased as k1 increased; Lc decreased as k1 increased while Lc increased as k2 increased; the principal bending moments and Lc increased as E increased. These calculations only analyzed the variations of the principal bending moments with a single parameter. To better characterize backfill mining, multiple parameters were analyzed to determine the effects of the principal bending moments of the main roof and Lc.(1)
Variations of the principal moments andLc with k1, k2, EFigure18 shows the variations of the principal bending moments and Lc with k1, k2, and E all taken into consideration. The bending moments of the main roof and Lc did not change when k1, k2 and E changed but the ratio of the three (k1/k2/E) remained constant, and the fracture position was unchanged.(2)
Changes in the maximum subsidenceω0 with various k1, k2, EFigure19 shows the changes in the maximum subsidence ω0 with various k1, k2 and E. The maximum subsidence ω0 of the main roof decreased gradually and the reduction extent significantly decreased when k1, k2 and E increased but with a constant ratio (k1/k2/E). The principal moments and fracture position of the main roof remained constant with the constant ratio (k1/k2/E), although the maximum subsidence ω0 changed significantly.Figure 18
Variation of the principal bending moments withk1/k2/E.Figure 19
Variation of the maximum subsidenceω0 with various k1, k2 and E.


### 4.8. Comparison of Model Conclusions

Under the condition of backfill mining, the mechanical model of thin plate was the most suitable to study the actual characteristics of overburden movement, while the rock beam model could only study the local mechanical characteristics in the mining area. Therefore, the model in this paper mainly studied the mechanical model of plate structure under the condition of backfill rather than the mechanical model of rock beam with defects.In general, the thin plate mechanical model with fixed support boundary around the overburden in the backfill area was widely used [31], so the comparison object selected in this paper was the thin plate mechanical model with fixed support boundary of four sides.The main differences were as follows:(1)
The elastic foundation coefficientk2 of backfill body, the elastic foundation coefficient k1 of the coal seam, the thickness h, and Young’s modulus E of main roof could affect the values and positions of the main roof principal moments, the distance of principal moments advancing coal wall. However, the traditional thin plate model with fixed support boundary could not get these useful conclusions.(2)
In this paper, it was found that the fracture location of the main roof was deep into the coal body area, while the traditional model showed that the first fracture location of the main roof was on the long side of the mining area and along the coal wall, which was independent of the boundary conditions, backfill body parameters, thickness, and Young’s modulus of the main roof.(3)
In particular, the weight relationship of the influencing factors of main roof fracture law was obtained; that is, when the ratio (k1/k2/E) was constant, the magnitude of the principal bending moments was constant, but the maximum subsidence changed significantly, which could not be obtained by the traditional plate structure model with fixed boundary.It could be seen that the conclusion of the traditional model had great limitations. The model and the conclusions of this paper were more practical and helpful to guide the engineering practice, which effectively made up for the shortcomings of the traditional model.


## 4.1. Fracture Criterion

(12)M1,M3=Mx+My2±Mx−My22+Mxy2,where M1 is the maximum principal bending moment; M3 is the minimum principal bending moment; Mx is the moment of x-component; My is the moment of y-component; and Mxy is the nodal twisting moment.(13)M1ijM3ij=Mxij+Myij2±Mxij−Myij22+Mxyij2,where (M1)ij is the nodal maximum principal moment; (M3)ij is the nodal minimum principal moment.Because the tensile strength of the rock is much less than the compressive strength, the principal moment of rock is compared with the ultimate bending moment of rock to judge whether the rock has failed or not. The maximum principal momentM1 and minimum principal moment M3 are obtained by equation (12), and the difference equation of equation (12) is equation (13). The moments (Mx)ij, (My)ij and (Mxy)ij in equation (13) can be obtained by equation (11). The extreme points of the principal moment are extracted from the cloud chart of the principal moment, which are shown in Figure 5.Figure 5
Principal moment feature whenk2 = 20 MN/m3. (a) Three-dimensional cloud chart of principal moment M1; (b) plane cloud chart of principal moment M1; (c) cloud chart of principal moment M3; (d) plane cloud chart of principal moment M3.
(a)
(b)
(c)
(d)


## 4.2. Effect of Elastic Foundation Coefficientk2 of Backfill Body in Backfill Mining

### 4.2.1. Basic Principles and Distribution Pattern of the Principal Moment

The elastic foundation coefficient of the backfill body is the key material factor to decide the backfill effect. Here, the backfill rate was 100%. The elastic foundation coefficientk2 was ∼20 MN/m3–220 MN/m3. When k2 = 20 MN/m3 and k2 = 180 MN/m3, cloud charts of the principal moment in backfill mining were plotted, as shown in Figures 5 and 6.Figure 6
Principal moment feature whenk2 = 180 MN/m3. (a) Three-dimensional cloud chart of principal moment M1; (b) plane cloud chart of principal moment M1; (c) three-dimensional cloud chart of principal moment M3; (d) plane cloud chart of principal moment M3.
(a)
(b)
(c)
(d)From Figures5 and 6, we can obtain the distribution of the principal moment during backfill mining: principal moments M1 and M3 in the middle area of the main roof directly overlying the backfill body were all positive bending moments; namely, the lower surface of the main roof bores the tensile stress and the upper surface of the main roof bores the compressive stress in this area. The nodal principal moments M3 on the long side and short side area with a certain distance to the coal wall (Lc) were negative bending moments; namely, the upper surface of the main roof bores the tensile stress and the lower surface of main roof bores the compressive stress in this area.According to the distribution characteristics of the principal moment in Figures5 and 6, the coordinates of the extreme points of the main roof’s principal moment were plotted, as shown in Figure 7.Figure 7
The coordinate of the extreme points of the main roof principal moment during backfill mining: (a)k2 = 20 MN/m3; (b) k2 = 180 MN/m3.
(a)
(b)When the value ofk2 was smaller (k2 = 20 MN/m3), and according to the distribution characteristics of the principal moment in Figure 5, the coordinates of the extreme points of the principal moment were plotted, as shown in Figure 7(a). The maximum absolute value of the principal moment of main roof in backfill area was at the midpoint with a coordinate of (0, 0), which was set to Mz, and Mz = M1|(0, 0); the maximum absolute value of the principal moment on the long side area with a certain distance to the coal wall (Lc) was the negative value of the minimum principal moment M3, which was set to Mc, the coordinates of which were (0, b + Lc) and (0, −b−Lc). The distance of the Mc advancing coal wall was set to Lc, and Mc=−M30,b+Lc=−M30,−b−Lc; the maximum absolute value of the principal moment on the short side area with a certain distance to the coal wall (Lc) was the negative value of the minimum principal moment M3, which was set to Md, and the coordinates of Md are (a + Ld, 0) and (−a−Ld, 0). The distance of Md advancing coal wall was set to Ld, and Md=−M3a+Ld,0=−M3−a−Ld,0.The coordinates of the extreme points of the principal moment from Figure6 were plotted when k2 was larger (k2 = 180 MN/m3), as shown in Figure 7(b). The maximum principal moment of main roof in backfill area was set to Mzc whose coordinates were (−a + Ldc, b−Lcc) and because of the symmetry, there were four points whose moments were all taken as MzcMzc−a+Ldc,b−Lcc>Mz=M10,0. In Figure 7(b), the maximum absolute value of the principal moment of the long side area with a certain distance to the coal wall (Lc) was set to Mc, and Mc=−M30,b+Lc=−M30,−b−Lc; the maximum absolute value of the principal moment of the short side area with a certain distance to the coal wall (Lc) wall was set to Md (Md=−M3a+Ld,0=−M3−b−Lc,0).As shown in Figures5 and 6, when the value of k2 became large (k2 = 180 MN/m3), the maximum principal moment Mz in the midpoint of the main roof transferred to the surrounding area and the points of the maximum principal moment were in the four corners of the study area. The maximum absolute value of the principal moment of main roof in the backfill area is denoted as Mzz, and Mzz meets the following relationship: Mzz=maxM10,0,Mzc−a+Ldc,b−LccBased on the above analysis, there appear to be three types of extreme points of the principal moments: the maximum absolute value of the principal momentMzz of the main roof in backfill area, the maximum absolute value of the principal moment Mc in the long side area, and the maximum absolute value of the principal moment Md in the short side area. Based on the above analysis and rock tensile properties, the lower surface of the middle area of main roof and the upper surface of long side and short side areas in front of the coal wall are the most susceptible to fracturing. Using the plate model on an elastic foundation with an elastic foundation boundary during backfill mining, the principles governing the fracturing of main roof are determined by the values of Mc, Mzz and Md. The main roof can only fracture when the principal moment is larger than Ms (the ultimate bending moment of the main roof).During backfill mining, the principal moments and the maximum subsidenceω0 of the main roof are related to k2, k1, h, E and μ. The principal movements and Lc were then calculated with variations in the above mentioned variables and analyzed in order to study the variations in the maximum subsidence of the main roof.


### 4.2.2. Effects of the Elastic Foundation Coefficientk2 of the Backfill Body

The following calculation examines the effect of the elastic foundation coefficientk2 on the principal bending moment and maximum subsidence ω0 of the main roof.(1)
The variation of the principal bending moment of the main roof andLc with the elastic foundation coefficient k2 of the backfill bodyFigure8 shows the variations of the principal bending moments of the main roof and Lc with the elastic foundation coefficient k2 of the backfill body. The backfill body significantly reduced the value of the maximum principal bending moment of the main roof and improved the stability of the main roof. The principal bending moments (Mc, Md and Mzz) continuously decreased with the increase in values of k2. The value of the principal bending moment Mzz was highly sensitive to k2, and comparatively the principal bending moment Md was less sensitive. The principal bending moments Mc and Md were finally almost equal with an increase in k2. In this situation, the long side and short side of the main roof almost fractured simultaneously when the maximum principal bending moment of the main roof increased to the ultimate bending moment Ms. There should be very little decrease in the principal bending moments when k2 increases to some extent.The initial fracture position of the main roof was in the long side area, which is located in front of the coal wall with aMc value greater than Ms. Lc continuously increased with the increase in k2. Specifically, when Mc was greater than Ms, the larger the values of k2, the further the distance between the fracture line of main roof and the coal wall (refer to Figure 8).(2)
Variation of the maximum subsidenceω0 of the main roof with an elastic foundation coefficient k2 of the backfill bodyFigure9 shows the variation of the maximum subsidence ω0 of the main roof with an elastic foundation coefficient k2 of the backfill body. The maximum subsidence ω0 of the main roof appeared to decrease with the increase in k2. The maximum subsidence should decrease very little when k2 increases to some extent.Using the above calculations and formulas, further study should be conducted to further describe the principal bending moments and the maximum subsidenceω0 of main roof with an elastic foundation coefficient k1 by varying the coal seam, thickness h of the main roof, Young’s modulus E of the main roof, and Poisson’s ratio μ of the main roof, in order to better and widely study fracture principles and subsidence of the main roof during backfill mining.Figure 8
Change in the principal bending moments andLc with elastic foundation coefficient k2 of backfill body.Figure 9
Change in the maximum subsidenceω0 of the main roof with variable k2.


## 4.2.1. Basic Principles and Distribution Pattern of the Principal Moment

The elastic foundation coefficient of the backfill body is the key material factor to decide the backfill effect. Here, the backfill rate was 100%. The elastic foundation coefficientk2 was ∼20 MN/m3–220 MN/m3. When k2 = 20 MN/m3 and k2 = 180 MN/m3, cloud charts of the principal moment in backfill mining were plotted, as shown in Figures 5 and 6.Figure 6
Principal moment feature whenk2 = 180 MN/m3. (a) Three-dimensional cloud chart of principal moment M1; (b) plane cloud chart of principal moment M1; (c) three-dimensional cloud chart of principal moment M3; (d) plane cloud chart of principal moment M3.
(a)
(b)
(c)
(d)From Figures5 and 6, we can obtain the distribution of the principal moment during backfill mining: principal moments M1 and M3 in the middle area of the main roof directly overlying the backfill body were all positive bending moments; namely, the lower surface of the main roof bores the tensile stress and the upper surface of the main roof bores the compressive stress in this area. The nodal principal moments M3 on the long side and short side area with a certain distance to the coal wall (Lc) were negative bending moments; namely, the upper surface of the main roof bores the tensile stress and the lower surface of main roof bores the compressive stress in this area.According to the distribution characteristics of the principal moment in Figures5 and 6, the coordinates of the extreme points of the main roof’s principal moment were plotted, as shown in Figure 7.Figure 7
The coordinate of the extreme points of the main roof principal moment during backfill mining: (a)k2 = 20 MN/m3; (b) k2 = 180 MN/m3.
(a)
(b)When the value ofk2 was smaller (k2 = 20 MN/m3), and according to the distribution characteristics of the principal moment in Figure 5, the coordinates of the extreme points of the principal moment were plotted, as shown in Figure 7(a). The maximum absolute value of the principal moment of main roof in backfill area was at the midpoint with a coordinate of (0, 0), which was set to Mz, and Mz = M1|(0, 0); the maximum absolute value of the principal moment on the long side area with a certain distance to the coal wall (Lc) was the negative value of the minimum principal moment M3, which was set to Mc, the coordinates of which were (0, b + Lc) and (0, −b−Lc). The distance of the Mc advancing coal wall was set to Lc, and Mc=−M30,b+Lc=−M30,−b−Lc; the maximum absolute value of the principal moment on the short side area with a certain distance to the coal wall (Lc) was the negative value of the minimum principal moment M3, which was set to Md, and the coordinates of Md are (a + Ld, 0) and (−a−Ld, 0). The distance of Md advancing coal wall was set to Ld, and Md=−M3a+Ld,0=−M3−a−Ld,0.The coordinates of the extreme points of the principal moment from Figure6 were plotted when k2 was larger (k2 = 180 MN/m3), as shown in Figure 7(b). The maximum principal moment of main roof in backfill area was set to Mzc whose coordinates were (−a + Ldc, b−Lcc) and because of the symmetry, there were four points whose moments were all taken as MzcMzc−a+Ldc,b−Lcc>Mz=M10,0. In Figure 7(b), the maximum absolute value of the principal moment of the long side area with a certain distance to the coal wall (Lc) was set to Mc, and Mc=−M30,b+Lc=−M30,−b−Lc; the maximum absolute value of the principal moment of the short side area with a certain distance to the coal wall (Lc) wall was set to Md (Md=−M3a+Ld,0=−M3−b−Lc,0).As shown in Figures5 and 6, when the value of k2 became large (k2 = 180 MN/m3), the maximum principal moment Mz in the midpoint of the main roof transferred to the surrounding area and the points of the maximum principal moment were in the four corners of the study area. The maximum absolute value of the principal moment of main roof in the backfill area is denoted as Mzz, and Mzz meets the following relationship: Mzz=maxM10,0,Mzc−a+Ldc,b−LccBased on the above analysis, there appear to be three types of extreme points of the principal moments: the maximum absolute value of the principal momentMzz of the main roof in backfill area, the maximum absolute value of the principal moment Mc in the long side area, and the maximum absolute value of the principal moment Md in the short side area. Based on the above analysis and rock tensile properties, the lower surface of the middle area of main roof and the upper surface of long side and short side areas in front of the coal wall are the most susceptible to fracturing. Using the plate model on an elastic foundation with an elastic foundation boundary during backfill mining, the principles governing the fracturing of main roof are determined by the values of Mc, Mzz and Md. The main roof can only fracture when the principal moment is larger than Ms (the ultimate bending moment of the main roof).During backfill mining, the principal moments and the maximum subsidenceω0 of the main roof are related to k2, k1, h, E and μ. The principal movements and Lc were then calculated with variations in the above mentioned variables and analyzed in order to study the variations in the maximum subsidence of the main roof.


## 4.2.2. Effects of the Elastic Foundation Coefficientk2 of the Backfill Body

The following calculation examines the effect of the elastic foundation coefficientk2 on the principal bending moment and maximum subsidence ω0 of the main roof.(1)
The variation of the principal bending moment of the main roof andLc with the elastic foundation coefficient k2 of the backfill bodyFigure8 shows the variations of the principal bending moments of the main roof and Lc with the elastic foundation coefficient k2 of the backfill body. The backfill body significantly reduced the value of the maximum principal bending moment of the main roof and improved the stability of the main roof. The principal bending moments (Mc, Md and Mzz) continuously decreased with the increase in values of k2. The value of the principal bending moment Mzz was highly sensitive to k2, and comparatively the principal bending moment Md was less sensitive. The principal bending moments Mc and Md were finally almost equal with an increase in k2. In this situation, the long side and short side of the main roof almost fractured simultaneously when the maximum principal bending moment of the main roof increased to the ultimate bending moment Ms. There should be very little decrease in the principal bending moments when k2 increases to some extent.The initial fracture position of the main roof was in the long side area, which is located in front of the coal wall with aMc value greater than Ms. Lc continuously increased with the increase in k2. Specifically, when Mc was greater than Ms, the larger the values of k2, the further the distance between the fracture line of main roof and the coal wall (refer to Figure 8).(2)
Variation of the maximum subsidenceω0 of the main roof with an elastic foundation coefficient k2 of the backfill bodyFigure9 shows the variation of the maximum subsidence ω0 of the main roof with an elastic foundation coefficient k2 of the backfill body. The maximum subsidence ω0 of the main roof appeared to decrease with the increase in k2. The maximum subsidence should decrease very little when k2 increases to some extent.Using the above calculations and formulas, further study should be conducted to further describe the principal bending moments and the maximum subsidenceω0 of main roof with an elastic foundation coefficient k1 by varying the coal seam, thickness h of the main roof, Young’s modulus E of the main roof, and Poisson’s ratio μ of the main roof, in order to better and widely study fracture principles and subsidence of the main roof during backfill mining.Figure 8
Change in the principal bending moments andLc with elastic foundation coefficient k2 of backfill body.Figure 9
Change in the maximum subsidenceω0 of the main roof with variable k2.


## 4.3. Effect of the Elastic Foundation Coefficientk1 of the Coal Seam in Backfill Mining

The elastic foundation coefficientk1 of the coal seam mainly reflects the limitation of the displacement of the main roof. We then analyzed variations of the principal bending moments and the maximum subsidence ω0 with an elastic foundation coefficient k1.(1)
Variations in the principal bending moments andLc with k1.Variations in the principal bending moments of the main roof and the distance to the coal wallLc with an elastic foundation coefficient k1 are presented in Figure 10. The principal bending moments (Mc, Md and Mzz) apparently increased with the increase in k1, and comparatively the principal bending moment Mzz had a smaller growth rate. When the principal moment Mc was greater than the ultimate moment Ms, the larger the value of k1 (the harder the coal seam), the closer the distance between the fracture line of the main roof and the coal wall (refer to Figure 10). Therefore, the influence of the boundary condition on the principal bending moment of main roof and Lc cannot be neglected. In other words, the influence of the boundary condition on the fracture position of main roof cannot be neglected.(2)
Variations in the maximum subsidenceω0 of the main roof with an elastic foundation coefficient k1 of coal seam.The change in the maximum subsidenceω0 with an elastic foundation coefficient k1 of the coal seam is presented in Figure 11. The maximum subsidence continuously decreased with an increase in k1, whereas the sensitivity of the maximum subsidence ω0 decreased with an increase in k1. The value of the maximum subsidence was reduced by approximately 10% with an increase in k1 (varying from 1000 to 8000 MN/m3).Figure 10
Change in the principal bending moments andLc with k1.Figure 11
Change in the maximum subsidenceω0 with k1.


## 4.4. Effect of the Thicknessh of the Main Roof

The thickness of the main roof can vary from mine to mine, depending on various circumstances (geology, previous mining, etc.). The effects of various thicknessesh were analyzed on the maximum subsidence ω0 and the principal bending moments in order to encompass a wide range of scenarios.(1)
Change in the principal bending moments andLc with hFigure12 graphically presents the variations of each of the three principal bending moments and Lc with various thicknesses of the main roof h (4–12 m). The principal bending moments (Mc, Md and Mzz) continuously increased with the increase in the roof thickness. Principal bending moment Mzz was highly sensitive to h, and comparatively the principal bending moment Md was less sensitive. The principal bending moments Mc and Md were almost equal when the thickness was relatively small (4–8 m). In such situation, the long side and short side area of the main roof almost fractured concurrently when the maximum principal bending moment of main roof increased to the ultimate bending moment Ms. Lc increased linearly with the increase in the roof thickness. Actually, the thicker the roof, the further the distance between the fracture line of main roof and the coal wall (refer to Figure 12).(2)
Changes in the maximum subsidence of the main roof with various thicknesseshFigure13 shows the change in the maximum subsidence ω0 of the main roof with various thicknesses (4–12 m). The maximum subsidence of the main roof continuously decreased with an increase in roof thickness, and the sensitivity of the maximum subsidence increased with an increase in the roof thickness.Figure 12
Change in the principal bending moments andLc with h.Figure 13
The maximum subsidenceω0 of main roof with roof thickness h.


## 4.5. Variations of Young’s ModulusE of the Main Roof in Backfill Mining

Young’s modulusE of the main roof was then varied to analyze the effects on the principal bending moment and maximum subsidence ω0.(1)
Changes in the principal bending moments andLc with various EFigure14 shows the changes in the principal bending moments and Lc with various Young’s modulus values E (10–75 GPa). The principal bending moments Mc and Md were almost equal when E was relatively small (10 GPa). In this condition, the long side and short side area of the main roof almost fractured concurrently when the maximum principal bending moment of the main roof increased to Ms. When the value of Young’s modulus reached 70 GPa, Mc and Mzz were nearly equal. In such situation, the long side and the middle area of the main roof almost fractured simultaneously when the principal bending moment was greater than Ms. The larger the E, the farther the Lc, and the farther the distance between the fracture line of main roof and the coal wall (refer to Figure 14). Based on the effect of Young’s modulus on the main bending moment, we can conclude that the greater Young’s modulus E value and the tensile strength must also be large within the main roof to ensure the main roof does not fracture.(2)
Change in maximum subsidenceω0 with various Young’s modulus E values of the main roofFigure15 shows the change in the maximum subsidence ω0 with an increase in Young’s modulus E of the main roof. Maximum subsidence of the main roof continuously decreased with an increase in E, and the sensitivity of ω0 decreased with an increase in E.Figure 14
Change in the principal bending moments andLc with E.Figure 15
The change in the maximum subsidenceω0 with various E.


## 4.6. Effects of Poisson’s Ratioμ on the Main Roof in Backfill Mining

Variations of the principal moments and maximum subsidenceω0 with Poisson’s ratio μ of main roof are discussed in the following calculations. Generally, the value of Poisson’s ratio varies from 0.1 to 0.3.(1)
Changes in the principal moments andLc with various Poisson’s ratio μ valuesFigure16 shows variations in principal moments and Lc with changes in Poisson’s ratio μ (0.10–0.30). The principal moments (Mc, Md and Mzz) gradually increased with the increase in μ, and Lc remained almost unchanged. Poisson’s ratio μ had little effect on principal moments and Lc; therefore, Poisson’s ratio would have little effect on the position of fracturing in the main roof.(2)
Variations of maximum subsidenceω0 of the main roof with μFigure17 shows how the maximum subsidence ω0 of the main roof varies with μ (0.10–0.30). The maximum subsidence slowly decreased with an increase in Poisson’s ratio. The maximum subsidence values ω0 when μ = 0.1 and μ = 0.3 were almost equal. Overall, Poisson’s ratio has little effect on the maximum subsidence of the main roof.Figure 16
The variation of principal bending moments andLc with μ.Figure 17
Variations of the maximum subsidenceω0 with μ.


## 4.7. Effects ofk1, k2, E in Backfill Mining

As observed through the above calculations, both the distribution of the principal bending moments of the main roof and the distance of the advancing coal seamLc were influenced by k1, k2 and E. The value of the principal bending moments decreased as k2 increased while the principal bending moments increased as k1 increased; Lc decreased as k1 increased while Lc increased as k2 increased; the principal bending moments and Lc increased as E increased. These calculations only analyzed the variations of the principal bending moments with a single parameter. To better characterize backfill mining, multiple parameters were analyzed to determine the effects of the principal bending moments of the main roof and Lc.(1)
Variations of the principal moments andLc with k1, k2, EFigure18 shows the variations of the principal bending moments and Lc with k1, k2, and E all taken into consideration. The bending moments of the main roof and Lc did not change when k1, k2 and E changed but the ratio of the three (k1/k2/E) remained constant, and the fracture position was unchanged.(2)
Changes in the maximum subsidenceω0 with various k1, k2, EFigure19 shows the changes in the maximum subsidence ω0 with various k1, k2 and E. The maximum subsidence ω0 of the main roof decreased gradually and the reduction extent significantly decreased when k1, k2 and E increased but with a constant ratio (k1/k2/E). The principal moments and fracture position of the main roof remained constant with the constant ratio (k1/k2/E), although the maximum subsidence ω0 changed significantly.Figure 18
Variation of the principal bending moments withk1/k2/E.Figure 19
Variation of the maximum subsidenceω0 with various k1, k2 and E.


## 4.8. Comparison of Model Conclusions

Under the condition of backfill mining, the mechanical model of thin plate was the most suitable to study the actual characteristics of overburden movement, while the rock beam model could only study the local mechanical characteristics in the mining area. Therefore, the model in this paper mainly studied the mechanical model of plate structure under the condition of backfill rather than the mechanical model of rock beam with defects.In general, the thin plate mechanical model with fixed support boundary around the overburden in the backfill area was widely used [31], so the comparison object selected in this paper was the thin plate mechanical model with fixed support boundary of four sides.The main differences were as follows:(1)
The elastic foundation coefficientk2 of backfill body, the elastic foundation coefficient k1 of the coal seam, the thickness h, and Young’s modulus E of main roof could affect the values and positions of the main roof principal moments, the distance of principal moments advancing coal wall. However, the traditional thin plate model with fixed support boundary could not get these useful conclusions.(2)
In this paper, it was found that the fracture location of the main roof was deep into the coal body area, while the traditional model showed that the first fracture location of the main roof was on the long side of the mining area and along the coal wall, which was independent of the boundary conditions, backfill body parameters, thickness, and Young’s modulus of the main roof.(3)
In particular, the weight relationship of the influencing factors of main roof fracture law was obtained; that is, when the ratio (k1/k2/E) was constant, the magnitude of the principal bending moments was constant, but the maximum subsidence changed significantly, which could not be obtained by the traditional plate structure model with fixed boundary.It could be seen that the conclusion of the traditional model had great limitations. The model and the conclusions of this paper were more practical and helpful to guide the engineering practice, which effectively made up for the shortcomings of the traditional model.


## 5. Engineering Practice

Changes in the principal moments of the main roof, the distance of theMc advancing coal wall Lc, and the maximum subsidence ω0 with k1, k2, E, h, μ and the ratio of k1/k2/E in backfill mining were examined to assist with forming the basis for engineering practice. To verify the theoretical calculations, the value of the elastic foundation coefficient k2 in the I26 backfill working face of Xingdong mine was determined. Subsidence of main roof and support pressure were then measured during mining.(1)
Calculating the elastic foundation coefficient of the backfill bodyThe tensile strength of the main roof was measured at 12.5 MPa and the ultimate bending moment of the main roof was 208 MN m in the I26 backfill working face. When the elastic foundation coefficient of the backfill body wask2 = 60 MN/m3, the maximum principal moment of main roof was Mc = 199 MN/m3, and the distance of the Mc advancing coal wall was Lc = 4 m (refer to Figure 8).To control the ground deformation and protect the main roof from fracturing, the value ofk2 of the backfill body should be about 56–60 MN/m3, while, using previous methods [31], the thin plate mechanical model with fixed support boundary, to calculate k2, the k2 would be 290 MN/m3. Clearly, the result calculated by the previous method was much larger than that calculated by the present method. From Figure 9 we know that the maximum subsidence of main roof given by theoretical calculation should be 44.5 mm when the value of k2 is 60 MN/m3.The I26 working face adopted the backfill bag of high-water material to backfill mining. In order to meet the production need of the working face and simultaneously ensure the convenience of backfilling the goaf, the size of the backfill bag was designed to be 15 × 2.2 × 4.8 m (length × width × height). The high-water material used in the I26 working face was a quick-setting high-water material. The material was prepared from the mixture of materials A and B, which solidify after mixing for 10–20 minutes. The elastic foundation coefficientk2 of the prepared high-water material was about 56–60 MN/m3.(2)
Field observationThe final main roof subsidence was about 48 mm, which was consistent with the value of theoretical calculations presented here.Mine pressure monitoring results showed that there was no periodic weighting during the mining of the I26 working face, which indicates that the main roof did not periodically fracture, and the borehole observation also indicated that the main roof did not fracture. Overall, it is feasible to determine the elastic foundation coefficientk2 of the backfilling material using the proposed theoretical calculation. The result calculated by the previous method was much larger than that calculated by the present method. Trying to achieve the larger elastic foundation coefficient of the backfilling material is much more costly and excessive for the given situation. Thus, the thin plate model on the elastic foundation with an elastic foundation boundary based on the actual boundary of the surrounding rock can better guide practice and can bring better economic benefits.


## 6. Conclusions

It is very important to establish a mechanical model that is in line with the actual characteristics of the project to effectively analyze and solve the mining engineering problems.In this paper, the thin plate model with elastic foundation boundary of overlying strata for backfill mining was established, and the calculation method was given. The results showed that the calculation method was feasible. Through this model, the fracturing law of the main roof was studied in detail, and new conclusions that the traditional model could not get were obtained, which made up the defects of the traditional model.(1)
The elastic foundation coefficientk1 of the coal seam, the elastic foundation coefficient k2 of backfill body, the thickness h, Young’s modulus E, and Poisson’s ratio μ of main roof could affect the values and positions of the main roof principal moments, the distance of principal moments advancing coal wall, and the maximum subsidence ω0; the influence of k2 on the principal moments and maximum subsidence ω0 of the main roof was the largest, while the influence of μ was the smallest; the values and positions of each of the principal moments remained unchanged but ω0 decreased significantly when k1, k2 and E increased and maintained the constant ratio (k1/k2/E).(2)
According to the variations in the principal moment, there were three types of initial fracturing of the main roof: fracturing in the long side ahead of the coal wall, fracturing in the long side ahead of the coal wall and the middle area of the main roof simultaneously, and fracturing in the long side ahead of the coal wall and short side ahead of the coal wall simultaneously, which was helpful to effectively monitor the fracturing position of the main roof in engineering practice.(3)
The elastic foundation coefficientk2 of the backfill body calculated by the thin plate model on the elastic foundation with elastic foundation boundary was far less than the value of k2 calculated by the thin plate model on an elastic foundation with fixed boundary. Using the I26 working face as verification of the calculation method, when the value of k2 was 56–60 MN/m3, the main roof did not undergo periodic weighting during mining. It provided a new and more reliable calculation method for determining the elastic foundation coefficient of the backfill body in order to protect the main roof from fracturing.(4)
The calculation model in this paper effectively made up for the defects and deficiencies of the traditional mechanical model, especially in the fracturing law of the main roof, fracturing conditions, and the weight relationship of the influencing factors. The conclusions promote the theoretical progress and had important reference significance for engineering practice.



---

*Source: 8906091-2020-03-17.xml*