Overburden and surface subsidence with slicing paste filling mining in thick coal seams

To overcome the difficulties of overburden failure and surface subsidence induced by the slicing mining of extra-thick coal seams, slicing filler paste is typically utilized. In this paper, a continuous curved beam mechanical model of paste filling mining was established by theoretical analysis against the background of the 3305 working face of Yangcheng Coal Mine, the main controlling factors of surface subsidence were analyzed using an orthogonal experiment method (OEM) and an analytic network process (ANP) coupled comprehensive assignment, and the vertical displacement distribution of the overlying rock under different filling intensity was simulated using numerical simulation software. The following main findings were obtained. First, the elastic modulus of the backfill is the primary regulating factor, as demonstrated by the results. Second, the greater the filling body’s elastic modulus, the more efficiently it carries the overburden load and minimizes the maximum surface subsidence. Third, the distance from the open-cut mine likewise reflects the degree of surface subsidence, with the amount of subsidence increasing as the distance from the mine increases.

Introduction

In recent years, coal resources have continued to dominate China’s energy structure, comprising more than 50% of the total share of all energy sources. However, the environmental issues caused by roof caving incidents and surface subsidence as a result of disturbances in coal mining remain severe. According to the relevant statistics (Wu and Li, 1995), it is typical in East China to mine coal under buildings, railways, and water bodies, which not only wastes resources but also compromises the integrity of the earth if the mining method is inappropriate. The backfill mining technique, which may prevent surface subsidence and maximize coal resources, is currently receiving widespread attention from domestic experts.

Numerous professionals and academics are now conducting research on roof overburden damage and surface subsidence, mostly through indoor experiments, theoretical analysis, numerical simulation, and field observations. Deng et al. (2015) suggested a novel mining technique of upward slicing filling mining for the geological circumstances of extra-thick coal seams using indoor experiments, investigated and researched this filling effect using similar material modeling, and conducted engineering practice experiments. Yang et al. (2015) utilized similarity theory to create comparable physical models to simulate the mining operation of underground ore deposits and monitor displacement changes. Liu et al. (2017a) constructed the theory of continuous curving beams for backfill mining using a large number of engineering applications and theoretical analysis as a starting point, revealing the mechanism for backfill control on the roof. Gao et al. (2017) established a defect-coupled constitutive model based on the theory of damage mechanics and showed that the piecewise nonlinear model more accurately depicts the effect of complex backfill defects on the stress–strain curve. Using geological data from Xiaotun Mine and numerical modeling, Zhang et al. (2011) assessed the influence of the filling ratio on surface subsidence. The results indicated that when the filling ratio increases, the control impact on surface subsidence improves. Kostecki and Spearing (2015) utilized finite-difference software to model the rheological characteristics of coal pillars in backfill mining and studied the link between the mechanical properties of the backfill material and the coal pillars under backfill mining settings. Numerous researchers have undertaken studies on the techniques and methods of backfilling, as well as the materials used in conjunction with the site, to avoid and manage roof-caused surface subsidence problems. Zhao et al. (2019) established an enhanced analytical method for evaluating the safety of vertically exposed backfill faces and maintaining the stability of the exposed backfill face. Thirukumaran et al. (2016) stated that the influence of the filling body on the surrounding rock is dependent on the mechanism, which consists of local support, surface support, and overall support. Gong et al. (2017) and Guo et al. (2018) investigated the physical parameters of filling with coal gangue as aggregate and the rationale and applicability of coal gangue mixture as backfill material.

During thick coal seam paste filling mining, the filling body in the mined-out region may handle the majority of the overburden load, which is one of the characteristics that separates it from caving. There is currently no consensus regarding important theories and application solutions for the site’s study status of the stress development law of the surrounding rock, overburdened rock, and surface movement deformation characteristics. The correlation between the strength of the filling material and the quantity of surface subsidence requires more investigation. To solve the problems of roof fall accidents and surface subsidence caused by thick coal seam mining, this paper analyzes the influencing factors of backfill, the effect of roof control, and predicts the amount of surface subsidence. Combining statistical principles and evaluation methods, the mining plan was finally chosen as using cement as the cementitious material, using tailing sand as the aggregate, mixed with fly ash and industrial slag as the filling material for backfill mining with the research background of the 3305 working face of Yangcheng Coal Mine.

Summary of the study area

Yangcheng Coal Mine is located in Jining City, Shandong Province, in the Liangshan coalfield in the southwestern block of the Luxi fault block and belongs to the North China sedimentary coalfield. Its stratigraphic distribution is Ordovician, Carboniferous, Permian, and Quaternary in order from the bottom to the top. There are two major coal-bearing strata in the study area, namely, the Upper Carboniferous Taiyuan Formation (C3t) and the Lower Permian Shanxi Formation (P1s). The 3# coal seams are the recoverable coal seams in the whole area, with the thickness ranging from 4 m to 9.5 m, and the average thickness being 7.5 m. The mine has an annual production capacity of 2.1 million t/a and a recoverable coal reserve of 133 million t in a near-horizontal seam. The 3305 working face is located on the west side of the third mining area, and the upper part is the 3303 working face goaf. The average coal seam thickness of the working face is 7.2 m, and the seam floor elevation is −870 m to −990 m. Due to the damage caused to the buildings and farmland of the nearby villages above the ground surface during the fully mechanized caving mining of the working face, the surface subsidence caused by mining was mitigated by using the paste filling mining method when extracting the coal seam of working face 3305; because the coal seam was thick, the whole coal seam could not be filled at one time, so the thick coal seam was extracted by using the slicing paste filling method. The working face layout and filling site are shown in Figure 1.

FIGURE 1

FIGURE 1. 3305 mining area working face layout and mining area filling site.

Mechanical modeling of roof overburden

Analysis of the mechanical model of continuous curved beams

During the mining operation, as the slurry is continuously filled into the goaf, the backfill material gradually solidifies and forms a stable support carrier, thereby limiting the sinking of the overlying rock layer, causing the rock layer to bend and sink, and preventing the irregular collapse of the roof (Shi et al., 2021). The pressure on the roof above the thick coal seam can be considered the uniform load (q₀), while the roof’s rocking beam is considered the elastic foundation beam. In addition, we consider the subterranean ceiling, the filling body, and the coal seam floor to be the elastic base. Figure 2 (Wang et al., 2020) depicts a structural mechanical model of thick coal slicing paste filling mining based on the aforementioned assumptions.

FIGURE 2

FIGURE 2. Mechanical model of slicing paste filling mining structure in thick coal seam.

Due to the supporting effect of the backfill body, the roof is more stable and less prone to fracture during the mining of thick coal seams with slicing paste filling. After mining damage, the roof is no longer supported on four sides and resembles a “simply supported beam” (Sandholtz et al., 2020). This structure’s mechanical model matches that of the upper and lower slicing mining surface. To assist further qualitative calculations, the load q is approximated as a uniform load q₀, and a schematic diagram of the force is presented in Figure 3.

FIGURE 3

FIGURE 3. Mechanical force sketch of the roof of the filling mining.

According to the Winkler foundation hypothesis, the pressure on any point on the foundation is positively correlated with the settlement at that point (Jena et al., 2020). Take the upper slicing face as an example, then,

$\left\{\begin{array}{cc}{p}_1\left(x\right)=k_1{w}_1\left(x\right)& \left(0\le x\le L\right)\\ p_2\left(x\right)=k_2{w}_2\left(x\right)& \left(x\le 0\right),\end{array}\right.(1)$

where p₁(x) is the support of the main roof by the foundation in the backfill area; p₂(x) is the support of the main roof by the foundation in the coal seam area; w₁(x) is the deflection of the main roof rock beam in the backfill area; w₂(x) is the deflection of the main roof rock beam in the coal seam area; k₁, k₂ is the elastic foundation coefficient determined by the elastic modulus; L is the strike length of the backfill area; and q₀ is the uniform load.

When mining the thick coal seam with slicing and filling, for the main roof deflection differential equation:

$\left\{\begin{array}{ll}EI\frac{{\mathrm{d}}^4{w}_1\left(x\right)}{d{x}^4}+k_1\left[w_1\left(x\right)-u_1\right]=q_0& \left(0\le x\le L\right)\\ EI\frac{{\mathrm{d}}^4{w}_2\left(x\right)}{d{x}^4}+k_2{w}_2\left(x\right)=q_0& \left(x\le 0\right),\end{array}\right.(2)$

where E is the roof elastic modulus; I is the moment of inertia of the roof section; and u₁ is the final subsidence amount of the immediate roof contacting the backfill body.

Set ${\beta }_1=\sqrt[4]{\frac{k_1}{4EI}},{\beta }_2=\sqrt[4]{\frac{k_2}{4EI}}$ and substitute the last formula:

$\left\{\begin{array}{ll}\frac{{\mathrm{d}}^4{w}_1\left(x\right)}{d{x}^4}+4{\beta }_1^4\left[w_1\left(x\right)-u_1\right]=\frac{q_0}{EI}& \left(0\le x\le L\right)\\ \frac{{\mathrm{d}}^4{w}_2\left(x\right)}{d{x}^4}+4{\beta }_2^4{w}_2\left(x\right)=\frac{q_0}{EI}& \left(x\le 0\right).\end{array}\right.(3)$

By solving the aforementioned equation we get

$\left\{\begin{array}{ll}{w}_1\left(x\right)=e^{-{\beta }_{1}x}\left(A\cos {\beta }_{1}x+B\sin {\beta }_{1}x\right)+\frac{q_0}{k_1}+u_1& \left(0\le x\le L\right)\\ w_2\left(x\right)=e^{{\beta }_{2}x}\left(C\cos {\beta }_{2}x+D\sin {\beta }_{2}x\right)+\frac{q_0}{k_2}& \left(x\le 0\right).\end{array}\right.(4)$

Since the mining method and the backfill body of the upper and the lower slicing face are the same, the amount of subsidence in the aforementioned Eq. 4 is the same. Combined with the boundary conditions, substituting x=0, $x=\pm \pi /2$ into Eq. 4, it can be determined that A, B, C, and D are as shown in Eq. 5:

$\left\{\begin{array}{l}A=-\left(\frac{q_0}{k_1}+u_1\right)\\ B=\frac{{\beta }_1-{\beta }_2}{{\beta }_1+{\beta }_2}\left(\frac{q_0}{k_1}+u_1\right)\\ C=\frac{{\beta }_1^2}{{\beta }_2^2}\left(\frac{q_0}{k_1}+u_1\right)\\ D=\frac{{\beta }_2-{\beta }_1}{{\beta }_1+{\beta }_2}\left(\frac{q_0}{k_1}+u_1\right).\end{array}\right.(5)$

Finally, the deflection curve equations of the upper slicing face and coal seam area are obtained as follows:

$\left\{\begin{array}{cc}{w}_1\left(x\right)=e^{-{\beta }_{1}x}\left(\frac{{\beta }_1-{\beta }_2}{{\beta }_1+{\beta }_2}\sin {\beta }_{1}x-\cos {\beta }_{1}x\right)\left(\frac{q_0}{k_1}+u_1\right)+\frac{q_0}{k_1}+u{}_1& \left(0\le x\le L\right)\\ w_2\left(x\right)=e^{{\beta }_{2}x}\left(\frac{{\beta }_2-{\beta }_1}{{\beta }_1+{\beta }_2}\sin {\beta }_{2}x+\cos {\beta }_{2}x\right)\left(\frac{q_0}{k_1}+u_1\right)\frac{{\beta }_1^2}{{\beta }_2^2}+\frac{q_0}{k_2}& \left(x\le 0\right).\end{array}\right.(6)$

In the same way, the deflection curve equation of the lower slicing face and the coal seam area can be obtained as follows (Jena et al., 2020):

$\left\{\begin{array}{cc}{w}_3\left(x\right)=e^{-{\beta }_{3}x}\left(\frac{{\beta }_3-{\beta }_4}{{\beta }_3+{\beta }_4}\sin {\beta }_{3}x-\cos {\beta }_{3}x\right)\left(\frac{q_0}{k_3}+u_1+u_2\right)+\frac{q_0}{k_3}+u{}_1+u_2& \left(0\le x\le L\right)\\ w_4\left(x\right)=e^{{\beta }_{4}x}\left(\frac{{\beta }_4-{\beta }_3}{{\beta }_3+{\beta }_4}\sin {\beta }_{4}x+\cos {\beta }_{4}x\right)\left(\frac{q_0}{k_3}+u_1+u_2\right)\frac{{\beta }_3^2}{{\beta }_4^2}+\frac{q_0}{k_4}& \left(x\le 0\right).\end{array}\right.(7)$

Main roof load and foundation coefficient

Main roof load

The main roof load during backfill mining can be regarded as the sum of the total pressure of the overlying strata, which is as follows:

$q={\gamma }_1{h}_1+{\gamma }_1{h}_1+\cdot \cdot \cdot +{\gamma }_m{h}_m.(8)$

Substituting the corresponding parameters of the overlying strata on the 3305 working face of Yangcheng Coal Mine 3# coal paste filling mining into the aforementioned formula, the main roof load q can be obtained as 22.5 MPa.

Foundation coefficient

Through the aforementioned analysis, it can be seen that the force and deformation of the roof are closely related to the foundation coefficient k (Chen and Park, 2020). When the coal seam is filled and mined, the elastic foundation is composed of the roof and the backfill body, which is regarded as an elastic body as a whole, and the constitutive relationship can be obtained as follows:

$\sigma =E\epsilon =E\frac{y}{h},(9)$

that is, the total compression of the immediate roof and the backfill body is as follows:

$\begin{array}{cc}{\sigma }_1=E_1\frac{y_1}{h_1},& {\sigma }_2=E_2\frac{y_2}{h_2}\end{array},(10)$

where: E₁ is the immediate roof elastic modulus; E₂ is the backfill body elastic modulus (Liu et al., 2017b); y₁ is the immediate roof compression; y₂ is the backfill body compression; h₁ is the immediate roof height (Mangal, 2021); and h₂ is the backfill body height.

According to the Winkler foundation assumption:

$\sigma =ky,(11)$

where $\sigma$ is the overburden load; k is the foundation coefficient; and y is the total compression of the roof and backfill body.

The foundation coefficient of the 3# coal seam lower slicing filling mining is similar to that of the upper slicing filling mining. Combining Eqs. 9–11, there is the following relationship between the foundation coefficient of the 0≤x≤L section and the foundation coefficient of the x≤0 section in the upper and lower slicing backfill mining:

$\begin{array}{c}\begin{array}{c}{k}_1=\frac{E_1{E}_2}{E_2{h}_1+E_1{h}_2}\\ k_2=\frac{E_1{E}_3}{E_3{h}_1+E_1{h}_3/2}\\ k_3=\frac{E_1{E}_2}{E_2{h}_1+E_1{h}_4}\\ k_4=\frac{E_1{E}_3}{E_3{h}_1+E_1{h}_3},\end{array}\end{array}(12)$

where k₁ is the foundation coefficient of upper slicing filling mining on 3# coal in the 0≤x≤L interval; k₂ is the foundation coefficient of upper slicing filling mining on 3# coal in the x≤0 interval; k₃ is the foundation coefficient of lower slicing filling mining on 3# coal in the 0≤x≤L interval; k₄ is the foundation coefficient of lower slicing filling mining on 3# coal in x≤0 interval; E₁ is the immediate roof elastic modulus; E₂ is the backfill body elastic modulus; E₃ is the coal seam elastic modulus; E₄ is the main roof elastic modulus; h₁ is the immediate roof height; h₂ is the height of the upper slicing backfill body; h₃ is the coal seam thickness; and h₄ is the backfill body height.

The relevant parameters of the upper face are as follows: E₁=6 GPa, h₁=4 m; E₂=0.5 GPa, h₂=3.24 m; E₃=3 GPa, h₃=7.2 m; E₄=21 GPa, and h₄=8 m (Ma et al., 2020). Substituting into the relevant parameters, we can obtain the foundation coefficient k₁=0.14 GN/m³, k₂=0.54 GN/m³.

The parameters of the immediate roof and the backfill body during the lower slicing mining are the same as those of the upper slicing mining. The overall height h₄ of the 3# coal backfill body is 6.48 m. Substituting the relevant parameters can obtain the foundation coefficient k₃=0.073 GN/m³, k₄=0.33 GN/m³.

Mechanical calculation of roof subsidence

Regarding the main roof subsidence characteristics of foundation in the filling area, from the aforementioned analysis, it can be seen that the main roof subsidence of the filling area of 3# coal upper and lower slicing paste filling mining is expressed as follows:

$\left\{\begin{array}{cc}\mathrm{T}\mathrm{h}\mathrm{e}\mathrm{u}\mathrm{p}\mathrm{p}\mathrm{e}\mathrm{r}:w_1\left(x\right)=e^{-{\beta }_{1}x}\left(\frac{{\beta }_1-{\beta }_2}{{\beta }_1+{\beta }_2}\sin {\beta }_{1}x-\cos {\beta }_{1}x\right)\left(\frac{q_0}{k_1}+u_1\right)+\frac{q_0}{k_1}+u{}_1& \left(0\le x\le l\right)\\ \mathrm{T}\mathrm{h}\mathrm{e}\mathrm{l}\mathrm{o}\mathrm{w}\mathrm{e}\mathrm{r}:w_3\left(x\right)=e^{-{\beta }_{3}x}\left(\frac{{\beta }_3-{\beta }_4}{{\beta }_3+{\beta }_4}\sin {\beta }_{3}x-\cos {\beta }_{3}x\right)\left(\frac{q_0}{k_3}+u_1+u_2\right)+\frac{q_0}{k_3}+u{}_1+u_2& \left(0\le x\le l\right).\end{array}\right.(13)$

The buried depth of the 3# coal seam is 900 m. The upper and the lower slicing mining areas have the same roof subsidence during the slicing paste filling mining; both are 0.36 m. We find the foundation coefficient k₁=0.14 GN/m³, k₂=0.57 GN/m³ in the upper slicing filling mining, and the foundation coefficient k₃=0.073 GN/m³, k₄=0.33 GN/m³, with the main roof load q=22.5 MPa in the lower slicing filling mining. The mechanical parameters of each influencing factor are shown in Table 1. We then substitute the relevant parameters in Table 1 into Eq. 13 and analyze the main roof subsidence’s characteristics of the upper and lower slicing mining filling areas (using MATLAB) and draw the fitting curve (using Origin) (Heinze et al., 2021). The schematic diagram of the influence curve of various factors on the main roof subsidence during upper slicing mining is shown in Figure 4, and the schematic diagram of the influence curve of various factors on the main roof subsidence during lower slicing mining is shown in Figure 5 (Deng et al., 2021).

TABLE 1

TABLE 1. Mechanical property parameters of different influencing factors.

FIGURE 4

FIGURE 4. Schematic diagram of the influence of various factors on the main roof subsidence (upper slicing mining). (A) Factor: main roof elastic modulus. (B) Factor: main roof thickness. (C) Factor: buried depth of the coal seam. (D) Factor: the subsidence of the roof contacting the backfill. (E) Factor: backfill elastic modulus. (F) Factor: coal elastic modulus.

FIGURE 5

FIGURE 5. Schematic diagram of the influence of various factors on the main roof subsidence (lower slicing mining). (A) Factor: main roof elastic modulus. (B) Factor: main roof thickness. (C) Factor: buried depth of the coal seam. (D) Factor: the subsidence of the roof contacting the backfill. (E) Factor: backfill elastic modulus. (F) Factor: coal elastic modulus.

In the upper and lower slicing mining, the influence of many sources on the principal roof subsidence has a clear, consistent tendency. As the distance from the coal wall (i.e., the distance of the open-off cut from the coal wall) grows, the backfill pressure first increases dramatically, peaks after a given distance, and then decreases gradually until it reaches the original rock stress (Yang et al., 2014). The following are the mechanical properties of the influencing factors’ parameters:

1) The modulus of elasticity and the thickness of the main roof has a certain influence on the pressure of the backfill body during coal seam filling and mining, which is manifested as follows: as the elasticity modulus and the thickness of the main roof increase, the maximum value of the pressure on the backfill body is increasingly further from the wall. Compared with the main roof elastic modulus, the influence of the main roof thickness is relatively significant.

2) The buried depth of the coal seam has a more obvious impact on the pressure of the backfill. The further the buried depth of the coal seam is from the surface, the greater the pressure of the backfill and the higher its peak value.

3) During coal seam filling and mining, the amount of subsidence when the roof touches the backfill also has a significant influence on the pressure on it. The greater the subsidence value, the higher the peak pressure of the backfill body.

4) During coal seam filling and mining, regarding the backfill body, the larger the elastic modulus, the larger the peak pressure on it, and as the peak pressure point gradually approaches the wall, the smaller the range of pressure on the backfill body.

5) Regarding the elastic modulus of the coal, as the elastic modulus of the coal gradually increases, the peak pressure of the backfill gradually decreases, but the pressure on the backfill is less affected.

Comprehensive risk assessment and analysis of main roof subsidence based on the adopted method

The alternative ranges of the main control factors are selected by the foundation coefficients and hydrogeological data, and the ranking of the importance of the weights is determined by using the comprehensive assignment method coupled with the analytic network process (ANP) and orthogonal experiment method (OEM) (Mosaad and Basheer, 2020). With the assessment methods based on expert experience and documentary data, the advantages of subjective analysis and objective analysis can be combined to improve prediction accuracy and credibility.

Analytic network process (ANP)

An overview of the analytic network process (ANP)

Thomas L. Saaty enhanced the analytic network technique based on the analytic hierarchy process (AHP) (Huang et al., 2020), which is an evaluation approach that combines qualitative and quantitative analysis. ANP can take into account the internal connections between nodes of various factor groups in a more thorough way than the classic AHP, which merely stresses the flaws of reciprocal impact across criteria levels (Kundu et al., 2021).

The weights are calculated as follows:

1) Construction of the judgment matrix. The control layer of the ANP has criterion elements B₁, B₂, B₃..., and the network layer under the control layer has element groups B_i1, B_i2, B_i3..., where i=1,2,....,n (Mirarabrazi and Navrodi, 2020). Take the control layer elements as the criteria and the network layer elements as the secondary criteria. The elements in the element group are compared according to their influence on B_ij, and the construction judgment matrix is obtained:

$W_{ij}=\left[\begin{array}{cccc}{w}_{\left(i1\right)}^{\left(j1\right)}& w_{\left(i1\right)}^{\left(j2\right)}& ...& w_{\left(i1\right)}^{\left(j{n}_j\right)}\\ w_{\left(i2\right)}^{\left(j1\right)}& w_{\left(i2\right)}^{\left(j2\right)}& ...& w_{\left(i2\right)}^{\left(j{n}_j\right)}\\ ...& ...& ...& ...\\ w_{\left(i{n}_i\right)}^{\left(j1\right)}& w_{\left(i{n}_i\right)}^{\left(j2\right)}& ...& w_{\left(i{n}_i\right)}^{\left(j{n}_j\right)}\end{array}\right].$

2) The ANP super cluster-weighted matrix and limit matrix (Camara et al., 2020). Combining the ranking vectors of the degree of mutual influence of the elements of the network layer, the supermatrix under the control element is obtained:

$W=\left[\begin{array}{cccc}{w}_{11}& w_{12}& ...& w_{1n}\\ w_{21}& w_{22}& ...& w_{2n}\\ ...& ...& ...& ...\\ w_{n1}& w_{n2}& ...& w_{nn}\end{array}\right].$

Each element of the aforementioned matrix is a non-negative matrix, and the columns of each supermatrix’s submatrices are normalized, but the entire matrix’s columns are not. The relevance of each element group is compared using the control level B as the criterion and the arbitrary element group B_i as the sub-criterion to generate the characteristic matrix:

$P=\left[\begin{array}{cccc}{p}_{11}& p_{12}& ...& p_{1n}\\ p_{21}& p_{22}& ...& p_{2n}\\ ...& ...& ...& ...\\ p_{n1}& p_{n2}& ...& p_{nn}\end{array}\right].$

Use the formula $\overline{W}=W·P$ to obtain the weighted supermatrix and use the formula $W^{\infty }=\underset{n\to \infty }{\lim }{\overline{W}}^n$ to obtain the limit matrix so as to obtain the importance ranking of each index.

Weight calculation by the analytic network process (ANP) model

According to the analysis of the mechanical model, we obtain six factors that affect the main roof subsidence: main roof elastic modulus; main roof thickness; buried depth of the coal seam; coal elastic modulus; backfill elastic modulus; and the subsidence of the roof contacting the backfill. Based on the attributes of the factors, they can be classified into three categorizations: roof structure; backfill characteristics; and coal lithology. A three-cluster ANP network model was developed as shown in Figure 6.

FIGURE 6

FIGURE 6. Three-cluster ANP network model.

Using group decision-making expert scoring, using 1–9 scales as the measurement criteria, and using (0,1) as the scoring range for the importance of the criteria or elements at the same level, pairwise comparison scores are obtained. The supermatrix composed of the judgment matrix is shown in Table 2, Table 3, Table 4, and Table 5.

TABLE 2

TABLE 2. Cluster-weight matrix.

TABLE 3

TABLE 3. Unweighted supermatrix.

TABLE 4

TABLE 4. Weighted matrix.

TABLE 5

TABLE 5. Limit matrix.

The maximum eigenroots and the corresponding eigenvectors are calculated according to the limit matrix, and then the eigenvectors are normalized to obtain the final weight of the roof subsidence amount evaluation index, as shown in Figure 7.

FIGURE 7

FIGURE 7. Distribution diagram of the weights for the main roof subsidence factors.

It can be seen from Figure 7 that the backfill elastic modulus and the subsidence of the roof contacting the backfill are the main controlling factors. They play a decisive role in the failure state of the overlying rock during sliced filling mining.

Verification experiment

Overview of the orthogonal experiment method (OEM)

Orthogonal experimental design and analysis methods are scientific calculation methods based on mathematical statistics, probability theory, and practical experience, using standardized orthogonal table design schemes and analyzing the results of precision and sensitivity. In the analysis of actual problems, the use of orthogonal tables to optimize the test plan can reduce the number of experiments while ensuring a certain level of accuracy, efficiently handling multi-factor optimization problems, and finding the main factors that have a significant impact on the test indicators among many factors.

The design process of the OEM (Xuan and Leung, 2011) is as follows: 1) determine the experimental factors and the number of levels; 2) select the applicable orthogonal table; 3) list the experimental plan and the experimental results; 4) perform range analysis and variance analysis on the experimental results; and 5) determine the optimal or better combination of factor levels.

Background of the orthogonal experiment method (OEM)

Weight calculation by the orthogonal experiment method (OEM)

Regarding at the six factors that affect the main roof subsidence in the 3305 coal face of Yangcheng Coal Mine, the primary and secondary degrees of each influencing factor are clarified by the OEM, and the orthogonal test table is designed as an example at 150 m behind the 3305 filling face of the lower slicing face. The main roof elastic modulus is selected as 15 GPa, 18 GPa, 21 GPa, 24 GPa, and 27 GPa; the main roof thickness is selected as 4 m, 6 m, 8 m, 10 m, and 12 m; the buried depth of the coal seam is selected as 700 m, 800 m, 900 m, 1000 m, and 1100 m; the subsidence of the roof contacting the backfill is selected as 0.05 m, 0.20 m, 0.35 m, 0.50 m, and 0.65 m; the backfill elastic modulus is selected as 0.1 GPa, 0.3 GPa, 0.5 GPa, 0.7 GPa, and 0.9 GPa; and the coal elastic modulus is selected as 1 GPa, 2 GPa, 3 GPa, 4 GPa, and 5 GPa. The orthogonal test table with six factors and five levels is established as shown in Table 6.

TABLE 6

TABLE 6. Orthogonal experiment table of various influencing factors.

We use the data in Table 6 to obtain the k_i value and the range value R of each factor through orthogonal calculation, and the results are shown in Table 7.

TABLE 7

TABLE 7. Orthogonal experiment table of various influencing factors.

The importance ranking was performed by finding the comprehensive average of the five levels (k₁, k₂, k₃, k₄, and k₅) of the six main control factors and conducting the range analysis to find the range value R corresponding to each factor. The order of importance positively reflects the importance of the factors, i.e., the first one is the most important and has the greatest weight, and so on. As can be seen from Table 7, the first two are the backfill elastic modulus and the subsidence of the roof contacting the backfill.

Determine the main control factors

The weight of contribution to the aim, as determined by the ANP expert scoring panel, is as follows: backfill elastic modulus > roof subsidence contacting the backfill > coal elastic modulus > buried depth of the coal seam > main roof thickness > main roof elastic modulus. The results were confirmed by the construction of orthogonal experiments based on objective data, and the ranking of the relevance of the main control factors derived by analysis of variance was compatible with the empirical analysis-based ANP group decision scoring approach. It can be determined that the elastic modulus of the backfill body is the main control factor for the subsidence of the backfill mining, and simulations were conducted for various orders of the backfill elastic modulus. In addition, the filling material that can reach the ideal modulus was chosen to handle the surface subsidence and destruction problem.

Numerical simulation of movement and deformation of overlying strata in the slicing paste filling mining of thick coal seams

Geomechanical parameters of the model

The numerical calculation model is established based on the engineering background of the geological mining conditions of the 3305 filling face of Yangcheng Coal Mine (Miranda et al., 2011). The simulated coal seam thickness is 7.2 m, the coal seam advancement distance is 750 m, and the working face width is 120 m. Due to the deep location of the coal seam, only the positions between 120 m above and below the coal seam (distance from the z-direction 240 m) are simulated, and the overburden gravity of 21.25 MPa is applied on the top. Considering the boundary effect factor, the simulated size is determined to be 1,500 m*350 m*463 m. The working face adopts the operation of aligning itself with the inclination, mining toward the direction, and filling while mining. The advancing amount is 1,000 m with 250-m coal pillars reserved at both ends, and 50-m protection coal pillars are reserved at both ends of the inclination. The simulation excavation is set with 20 excavation steps, where each step is 50 m, and a total of 109,224 grids are generated by the finite difference method of FLAC3D. From the 60 observation points, we selected five representative points (A, B, C, D, and E) with different orientations as sample data collection points. Taking the orientation of the model as the X axis, the inclination as the Y axis, and the height as the Z axis, an XOY space rectangular coordinate system is established. The coordinates of A, B, C, D, and E are as follows: A (300;100;400); B (500;150;360); C (700;200;320); D (900;280;280); and E (1,100;300;240).The geomechanical parameters of coal and rock are shown in Table 8, and the numerical model is shown in Figure 8.

TABLE 8

TABLE 8. Geomechanical parameters of coal and rock formations.

FIGURE 8

FIGURE 8. Numerical calculation model.

Overlying strata and surface deformation characteristics of thick coal seam slicing paste filling mining under different backfill strengths

Paste filling materials can be calculated by the Salamon formula double-yield plasticity model (Du et al., 2020; Wang and Du, 2020; Yang et al., 2020). From the conclusion of the previous section, it can be seen that the strength of the backfill is a key factor affecting the filling effect of the paste (Liu et al., 2017c). In order to verify the effect of backfill strength on overburden and surface control, numerical simulation was conducted to analyze the overburden and surface subsidence by comparing three cases of 0.3 GPa, 0.5 GPa, and 1.0 GPa backfill strength at a 95% fill rate.

We compared three different backfill strengths of 0.3 GPa, 0.5 GPa, and 1.0 GPa. We took the vertical displacement distribution of overburden strata after all the advancing of the 3# coal 3305 working face as an example for comparison, as shown in Figure 9, and the surface subsidence comparison graph is shown in Figure 10.

FIGURE 9

FIGURE 9. Vertical displacement distribution of overlying strata under different backfill strengths. (A) Cloud image of vertical displacement distribution of 1.0 GPa 3# coal overlying strata. (B) Cloud image of vertical displacement distribution of 0.5 GPa 3# coal overlying strata. (C) Cloud image of vertical displacement distribution of 0.3 GPa 3# coal overlying strata.

FIGURE 10

FIGURE 10. Diagram of the overburden movement process during slicing mining. (A) 1.0 GPa roof subsidence of the upper slicing face. (B) 1.0 GPa roof subsidence of the lower slicing face. (C) 0.5 GPa roof subsidence of the upper slicing face. (D) 0.5 GPa roof subsidence of the lower slicing face. (E) 0.3 GPa roof subsidence of the upper slicing face. (F) 0.3 GPa roof subsidence of the lower slicing face.

Five points with different orientations (A, B, C, D, and E) were selected from the 60 observation points in FLAC3D numerical simulation to observe the amount of roof subsidence corresponding to the filling of the slicing mining with different elastic moduli of 0.3 GPa, 0.5 GPa, and 1.0 GPa at a 95% filling rate. With the advancement of the 3305 working face, the overlying rock layer gradually sinks on the original foundation, and the amount of subsidence increases with the advancement of mining. The upper slicing face reaches its peak at an excavation step distance of about 800 m, and the lower slicing face reaches its peak at an excavation step distance of about 1,500 m; until the direction of mining is infinite, the amount of subsidence tends to be stable. We observed the change of the effect of surface subsidence and overburden damage at five points (A, B, C, D, and E) using backfill with different elastic moduli. The subsidence of the roof with mining damage was ranked 0.3 GPa>0.5 GPa>1.0 GPa. In the maximum settlement curves of the upper and lower strata, 1.0 GPa<0.5 GPa<0.3 GPa, and the peak of the main roof subsidence of the upper slicing face is greater than that of the lower slicing face.

Analysis of mine pressure monitoring data

According to the mine monitoring pressure data before and after backfill mining at the 3305 working face, a comparative analysis was conducted. With the advancement of mining, the pressure on the protective coal pillars showed different changes, and the mine monitoring pressure graphs before and after filling were obtained as shown in Figure 11.

FIGURE 11

FIGURE 11. Rock pressure monitoring on the 3305 working face. (A) Rock pressure monitoring before filling. (B) Rock pressure monitoring after filling.

From Figure 11A, it can be seen that in the unfilled area, the mine pressure is obvious and the pillar is under greater force. During the mining process, the mine pressure started to decrease at 10 m of excavation because the coal body was extracted one after another. When the excavation reaches 40 m, the periodic incoming pressure appears, and comparing this part of mine pressure with the previous numerical simulation, it can be found that it basically matches. According to Figure 11B, it can be seen that, after filling, the mine incoming pressure shows a large and stable decrease, and the increase tends to be stable in the later period. The peak and average values of the pillar pressure are significantly reduced, which reduces the risk of rock impact and roof collapse and deformation. Therefore, it can be concluded that slicing paste filling mining with cement, tailing sand, fly ash, and industrial slag as the filling materials for thick coal seams can effectively support the roof and reduce the amount of subsidence.

Conclusion

1) The mechanical model of the roof elastic foundation beam model was constructed using mechanical knowledge, and the deflection curve equations in the filling area and coal body area were established. Through the force analysis of the roof, the range of factors that affect the amount of roof sinking was obtained. The relationship curves between the roof elastic modulus, the roof thickness, the buried depth of the coal seam, the elastic modulus of backfill, the coal elastic modulus, and the subsidence were drawn.

2) From the expression of the main roof subsidence in the filling area of 3# coal upper and lower slicing paste filling mining and the geological data of the 3305 working face of Yangcheng Coal Mine, six influencing factors affecting the amount of main roof subsidence were inductively analyzed, and the weights were derived using the subjective group decision method of the ANP and verified by the OEM based on objective data, and the importance ranking coinciding with the weight ranking was derived. The results show that the backfill elastic modulus is the main controlling factor.

3) FLAC3D was used to model backfill elastic moduli of 0.3 GPa, 0.5 GPa, and 1 GPa during excavation while maintaining a constant filling rate. The simulation revealed the effect of the filling body’s elastic modulus on the overburden and surface sinking. The results demonstrate that increases the elastic modulus of the backfill reduces the degree of overburden damage and the quantity of surface subsidence, indicating that increasing the elastic modulus of the filling body can effectively enhance the effect while maintaining the filling rate.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary material.

Author contributions

All authors listed have made a substantial, direct, and intellectual contribution to the work and approved it for publication. QZ was responsible for writing the paper, deriving the model, and translating and embellishing it. CW directed the structure of the thesis. LP was responsible for the literature proofreading.

Funding

This research was funded by the National Emergency Management System Construction Project (grant 20VYJ061) and the Construction and Empirical Research on Early Warning Index System of Major Engineering Safety Risks Based on Optimal Control Theory, National Natural Science Foundation of China (grant 71271031).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Camara, M., Jamil, N. R., Abdullah, A. F. B., Hashim, R. b., and Aliyu, A. G. (2020). Economic and efficiency based optimisation of water quality monitoring network for land use impact assessment. Sci. Total Environ. 737, 139800. doi:10.1016/j.scitotenv.2020.139800

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, Y., and Park, Y. H. (2020). A Helmholtz resonator on elastic foundation for measurement of the elastic coefficient of human skin. J. Mech. Behav. Biomed. Mat. 101, 103417. doi:10.1016/j.jmbbm.2019.103417

PubMed Abstract | CrossRef Full Text | Google Scholar

Deng, X. J., Liu, H., Feng, J., Li, Y., de Wit, B., and Yang, Y. (2021). A continuous Conversion Technology from Top coal caving into multi-slice mining for thick Bifurcated coal seams. Geotech. Geol. Eng. (Dordr). 39, 5983–5992. doi:10.1007/s10706-021-01898-6

CrossRef Full Text | Google Scholar

Deng, X. J., Zhang, J. X., Huang, P., Zhang, Q., and Hao, X. J. (2015). Analysis of roof movement characteristics in upward slicing and filling mining in ultra-thick coal seams. J. Chin. Coal. Soc. 40, 994–1000.

Google Scholar

Du, F., Wang, K., Zhang, X., Xin, C., Shu, L., and Wang, G. (2020). Experimental study of coal–gas outburst: Insights from coal–rock structure, gas pressure and adsorptivity. Nat. Resour. Res. 29, 2481–2493. doi:10.1007/s11053-020-09621-7

CrossRef Full Text | Google Scholar

Gao, R., Zhou, K., and Yang, C. (2017). Damage mechanism of composite cemented backfill based on complex defects influence. Materwiss. Werksttech. 48, 893–904. doi:10.1002/mawe.201600722

CrossRef Full Text | Google Scholar

Gong, P., Ma, Z., Zhang, R. C., Ni, X. Y., Liu, F., and Huang, Z. (2017). Surrounding rock deformation mechanism and control Technology for Gob-side entry Retaining with fully mechanized gangue backfilling mining: A case study. Shock Vib. 3, 1–15. doi:10.1155/2017/6085941

CrossRef Full Text | Google Scholar

Guo, P., Meng, Z. G., Ning, X. Y., and Zhang, R. R. (2018). An experimental investigation on the mechanical properties of gangue Concrete as a Roadside support body material for backfilling Gob-side entry Retaining. Adv. Mater. Sci. Eng. 17, 1–11. doi:10.1155/2018/1326053

CrossRef Full Text | Google Scholar

Heinze, T., Frank, S., and Wohnlich, S. (2021). FSAT - A fracture surface analysis toolbox in MATLAB to compare 2D and 3D surface measures. Comput. Geotech. 132, 103997. doi:10.1016/j.compgeo.2020.103997

CrossRef Full Text | Google Scholar

Huang, Z., Le, T., Gao, Y., Yao, X., Wang, H., Zhao, W., et al. (2020). Safety assessment of Emergency Training for industrial accident Scenarios based on analytic Hierarchy process and Gray-Fuzzy comprehensive assessment. IEEE Access 8, 144767–144777. doi:10.1109/access.2020.3013671

CrossRef Full Text | Google Scholar

Jena, S. K., Chakraverty, M., and Ma, L. (2020). Vibration and buckling characteristics of nonlocal beam placed in a magnetic field embedded in Winkler–Pasternak elastic foundation using a new refined beam theory: An analytical approach. Eur. Phys. J. Plus 135, 164–218. doi:10.1140/epjp/s13360-020-00176-3

CrossRef Full Text | Google Scholar

Kostecki, T., and Spearing, A. (2015). Influence of backfill on coal pillar strength and floor bearing capacity in weak floor conditions in the Illinois Basin. Int. J. Rock Mech. Min. Sci. (1997). 76, 55–67. doi:10.1016/j.ijrmms.2014.11.011

CrossRef Full Text | Google Scholar

Kundu, A., Dutta, D., Patel, N. R., Denis, D. M., and Chattoraj, K. K. (2021). Evaluation of Socio-Economic Drought risk over Bundelkhand region of India using analytic Hierarchy process (AHP) and Geo-Spatial techniques. J. Indian Soc. Remote Sens. 49, 1365–1377. doi:10.1007/s12524-021-01306-9

CrossRef Full Text | Google Scholar

Liu, E. Y., Zhang, Q. l., Feng, Y., and Zhao, J. w. (2017). Experimental study of static and dynamic mechanical properties of double-deck backfill body. Environ. Earth Sci. 76, 689–6897. doi:10.1007/s12665-017-7011-0

CrossRef Full Text | Google Scholar

Liu, J. G., Zhao, J. W., and Yang, H. Z. (2017). Research on the time and space characteristics of elastic foundation beam in filling mining. Coal Sci. Technol. 45, 41–47.

Google Scholar

Liu, J., Sui, W., and Zhao, Q. (2017). Environmentally sustainable mining: A case study of intermittent cut-and-fill mining under sand aquifers. Environ. Earth Sci. 76, 562. doi:10.1007/s12665-017-6892-2

CrossRef Full Text | Google Scholar

Ma, Q., Tan, Y., Liu, X., Gu, Q., and Li, X. (2020). Effect of coal thicknesses on energy evolution characteristics of roof rock-coal-floor rock sandwich composite structure and its damage constitutive model. Compos. Part B Eng. 198, 108086. doi:10.1016/j.compositesb.2020.108086

CrossRef Full Text | Google Scholar

Mangal, A. (2021). A study of stability investigation of immediate roof for extraction of thick coal seam in India. J. Ins. Eng. (India) Ser. D., 102, 203–221.

CrossRef Full Text | Google Scholar

Miranda, T., Dias, D., Eclaircy-Caudron, S., Gomes Correia, A., and Costa, L. (2011). Back analysis of geomechanical parameters by optimisation of a 3D model of an underground structure. Tunn. Undergr. Space Technol. 26, 659–673. doi:10.1016/j.tust.2011.05.010

CrossRef Full Text | Google Scholar

Mirarabrazi, J., and Navrodi, I. H. (2020). Identifying optimal location of ecotourism sites by analytic network process and genetic algorithm (GA): (Kheyroud forest). IJEST 17, 2583–2592.

CrossRef Full Text | Google Scholar

Mosaad, S., and Basheer, A. A. (2020). Utilizing the Geophysical and hydrogeological data for the assessment of the Groundwater Occurrences in Gallaba plain, western Desert, Egypt. Pure Appl. Geophys. 177, 3361–3382. doi:10.1007/s00024-019-02414-x

CrossRef Full Text | Google Scholar

Sandholtz, S. H., Kannan, D., Beltran, B. G., and Spakowitz, A. J. (2020). Chromosome structural mechanics Dictates the local Spreading of Epigenetic Marks. Biophys. J. 119, 1630–1639. doi:10.1016/j.bpj.2020.08.039

PubMed Abstract | CrossRef Full Text | Google Scholar

Shi, X., Zhou, H., Sun, X., Cao, Z., and Zhao, Q. (2021). Floor damage mechanism with cemented paste backfill mining method. Arab. J. Geosci. 14 (2021), 80–89. doi:10.1007/s12517-020-06368-6

CrossRef Full Text | Google Scholar

Thirukumaran, S., Indraratna, B., Brown, & E. T., and Kaiser, P. K. (2016). Stability of a rock block in a Tunnel roof under constant normal Stiffness conditions. Rock Mech. Rock Eng. 49, 1587–1593. doi:10.1007/s00603-015-0770-6

CrossRef Full Text | Google Scholar

Wang, H., Wang, E., Li, Z., Shen, R., Li, B., Zhang, Q., et al. (2020). Study and application of dynamic inversion model of coal seam gas pressure with drilling. Fuel 280, 118653. doi:10.1016/j.fuel.2020.118653

CrossRef Full Text | Google Scholar

Wang, K., and Du, F. (2020). Coal-gas compound dynamic disasters in China: A review. Process Saf. Environ. Prot. 133, 1–17. doi:10.1016/j.psep.2019.10.006

CrossRef Full Text | Google Scholar

Wu, K., and Li, B. (1995). Energy development in China: National policies and regional strategies. Energy Policy 23, 167–178. doi:10.1016/0301-4215(95)91420-h

CrossRef Full Text | Google Scholar

Xuan, W., and Leung, D. (2011). Optimization of biodiesel production from camelina oil using orthogonal experiment. Appl. Energy 88, 3615–3624. doi:10.1016/j.apenergy.2011.04.041

CrossRef Full Text | Google Scholar

Yang, R. S., Zhu, Y., Li, Y., Li, W., and Lin, H. (2020). Coal pillar size design and surrounding rock control techniques in deep longwall entry. Arab. J. Geosci. 13, 453. doi:10.1007/s12517-020-05454-z

CrossRef Full Text | Google Scholar

Yang, W., Lin, B. Q., and Xu, J. T. (2014). Gas outburst affected by original rock stress direction. Nat. Hazards (Dordr). 72, 1063–1074. doi:10.1007/s11069-014-1049-z

CrossRef Full Text | Google Scholar

Yang, Z. Q., Zhai, S. H., Gao, Q., and Li, M. H. (2015). Stability analysis of large-scale stope using stage subsequent filling mining method in Sijiaying Iron Mine. J. Rock Mech. Geotech. Eng. 7, 87–94. doi:10.1016/j.jrmge.2014.11.003

CrossRef Full Text | Google Scholar

Zhang, Q. L., Zhou, H. Q., and Bo, J. B. (2011). Research and Practice on stability of overlying strata in paste filling mining. J. Min. Saf. Eng. 28, 279–282.

Google Scholar

Zhao, X., Fourie, A., and Qi, C. C. (2019). An analytical solution for evaluating the safety of an exposed face in a paste backfill stope incorporating the arching phenomenon. Int. J. Min. Metall. Mat. 26, 1206–1216. doi:10.1007/s12613-019-1885-7

CrossRef Full Text | Google Scholar