a mining subsidence prediction model under thick loose layer and its parameter inversion method

How to cite item Li, J., Yu, X., Chen, D., & Fang, X. (2021). Research on the establishment of a mining subsidence prediction model under thick loose layer and its parameter inversion method. Earth Sciences Research Journal, 25(2), 215-223. DOI: https://doi.org/10.15446/esrj.v25n2.79537 Most of the coal mining in China is underground, which will inevitably cause surface deformation and trigger a series of geological disasters. Therefore, it is essential to find a suitable method to forecast the ground sinking caused by underground mining. The most commonly used prediction model in China is the probability integral model (PIM). But when this model is used in the geological condition of mining under thick loose layers, the predicted edge of the sinking basin will converge faster than the actual measured sinking situation. A geometric model (GM) with a similar model shape as the PIM but with a larger boundary value was established in this paper to solve this problem. Then an improved cuckoo search algorithm (ICSA) was proposed in this paper to calculate the GM parameters. The stability and reliability of the ICSA were verified through a simulated working face. At last, the ICSA, in combination with the GM and the PIM, was used to fit 6 working faces with the geological mining condition of thick loose layers in the Huainan mining area. The results prove that GM can solve the above-mentioned PIM problem when it is used in geological mining conditions of thick loose layers. And it was obtained through comparative analysis that the GM and the PIM parameters can take the same value except for the main influence radius. ABSTRACT Research on the establishment of a mining subsidence prediction model under thick loose layer


Introduction
Although the continuous mining of the coal resources provides a strong guarantee for the rapid development of China's economy. However, on the one hand, with the continuous consumption of coal resources and the continuous acceleration of urbanization, the amount of the coal seams pressed under buildings, railways, and water is increasing. On the other hand, the surface subsidence and the massive structure damage caused by underground mining have been very serious in China (Xuan & Xu, 2014;Zhao et al., 2016;Li et al., 2019;Xu et al., 2019). Therefore, it is very important to predict and fully grasp the laws of the ground deformation and its deformation caused by underground mining, which can provide guidance for controlling the destruction of surface protection objects. However, the process of the surface subsidence is complicated and different geological mining conditions will cause different surface subsidence laws . For example, the mining subsidence laws under the geological condition of thick loose layers have special characteristics compared with thinner loose layers and no loose layer. Such as the maximum subsidence value may be greater than the normal thickness of the coal seam and it is apparently smaller of the predicted value than the measured data at the boundary of the subsidence basin Zhou et al., 2015).
Coal seams buried under thick loose layers are widely distributed in China. In recent years, most research on mining subsidence under the geological condition of thick domestic and foreign scholars has focused on the mechanism of rock formation movement. Zhu et al. (2020) confirmed that the transition layers of the bedrock layers and the loose layers has a support effect on mining subsidence. Liu et al. (2020) found that the thick loose layer formed a bending collapse zone, fracture development zone and caving zone during the mining process through simulation analysis, and derived a prediction model based on elastic mechanics. Xu et al. (2020) studied the movement laws of the overlying bedrock while mining shallow buried coal seams under thick loose layers. Yu et al. (2020) studied the rock-soil structure above the coal seam which is buried under thin bedrock layer and thick loose layer. Wang et al. (2019) analyzed the stability of the loose layers and the bedrock layers during mining. Ma et al. (2017) discussed the prediction methods for mining subsidence when coal pillars are solid backfilled under thick loose layers. Yang & Xia (2013) studied the ground deformation laws caused by mining under thick loose layers and thin bedrock layers by analyzing measured data. Dai et al. (2011) researched the influence of the thickness of the loose layers on mining subsidence through numerical simulation. Zhou et al. (2018) proposed a combined prediction model based on soil mechanics and random medium theory. Chang et al. proposed a mining subsidence prediction model used in subcritical mining under thick loose layers with considering the different subsidence and movement mechanisms of the bedrock layers and the thick loose layers.
But the model can only be suitable for subcritical mining under thick loose layers (Chang et al., 2015). The advantage of the mining subsidence prediction model which is based on the mechanism of rock formation can fundamentally find the law of rock formation. It is not easy to establish the model and calculate the model parameters in actual application. So they all have not been widely used. The most commonly used mining subsidence prediction model is always the PIM which is developed from the random medium theory Zheng et al., 2019;Diao et al., 2016;Fan et al., 2015). However, due to the limitation of the shape of the model curve, the problem of small boundary value is much more obvious while the PIM is applied to predict surface subsidence caused by mining under thick loose layers, and this problem cannot be solved by adjusting the model parameters. Therefore, it is necessary to establish a model with simple parameters and suitable for surface subsidence prediction caused by mining under thick loose layers. An accurate prediction model is very important to mining subsidence prediction, a set of accurate parameters to the model is also important (Sehn et al., 2015). At present, the mostly used parameters inversion methods in mining subsidence prediction area can be divided into direct optimization methods (Shen et al., 2015;Guo & Wang, 2000) and intelligent optimization algorithms (Zha et al., 2011;Wang et al., 2018). The most commonly used of direct optimization methods is the modular vector method (MVM). The MVM has the advantages of easy programming and is suitable for working faces of any shape. But it needs a set of initial values of the parameters while using, and use them as a basis to follow the valley line or the ridge line to accelerate to the best values. The MVM has a strong dependence on the initial value of parameters. If the initial value is not accurate, the calculation result will deviate from the reasonable range. Since the PIM is an official forecasting method in China, many mining faces have a set of PIM parameters. Moreover, the subsidence prediction parameters of the working face with similar geological conditions are similar, so the PIM parameters of the mined working face with similar geological conditions can be used as the initial value of the PIM parameter of the working face to be solved. But on the one hand, the initial values of the parameters obtained in this way may be unreliable. On the other hand, it is difficult to obtain the initial values of the parameters of new prediction models. Therefore, the use of the MVM was restricted. The intelligent optimization algorithm uses random global search within a certain range by simulating the behavioral characteristics of creatures such as foraging and reproduction instead of deterministic search of the modular vector method, so it can explore the global optimal value without initial value dependence. The most commonly used method of intelligent optimization algorithms is the genetic algorithm (GA) (Du et al., 2014;Xing et al., 2021;Li et al., 2017). However, the accuracy of the GA to calculate the parameters of the mining subsidence prediction model is not ideal.
At the beginning of this paper, the GM to predict the surface subsidence caused by limited coal seams is established. And it is found that by adjusting the parameters of the GM, the middle part of its model curve can be basically consistent with the PIM, while the boundary value predicted is larger than the PIM. Therefore, it is inferred that the shape of the GM curve is more in line with the shape of the main section of the subsidence basin which is caused by mining under thick loose layers. Then in order to calculate the GM parameters, an improved cuckoo search algorithm-ICSA-is proposed, and a simulated working surface is used to verify the stability and accuracy of the ICSA. Finally, 6 working faces buried under thick loose layers were used to compare the accuracy of the GM and the PIM. By comparing the fitting effects of the GM and the PIM in 6 working faces, it can be concluded that the GM is more suitable than the PIM in predicting surface subsidence caused by mining under thick loose layers. and the relationship between the GM parameters and the PIM parameters is drawn.

Semi-infinite mining subsidence prediction model
Take the open-off cut of the coal seam as the origin, the direction from the origin to the inclination direction of the goaf as the positive direction of the y-axis and the direction from the origin to the strike of the goaf as the positive direction of the x-axis to establish the calculation coordinate system. Then according to the PIM, the sinking value of the point with abscissa x on the strike main section of the subsidence basin caused by a semi-infinite mining coal seam with the thickness of m and the dip angle of α can be calculated by Equation 1.
Where W 0 represents the max subsidence value and can be calculated by W m cos , q PIM is a parameter in the PIM named subsidence coefficient whose value is mainly determined by the lithology of the overlying rock, and its value is generally less than 1. r PIM is another parameter in the PIM named main influence radius.
It can be seen from Equation 2 that the main section subsidence curve caused by a semi-infinite mining coal seam is in the shape of "S". And the .
, q GM is the subsidence coefficient of the GM. r GM is the main influence radius of the GM.

Limited mining subsidence prediction model
The sinking value of the point with abscissa x on the strike main section of the subsidence basin caused by a limited mining coal seam with the length of D 3 can be regarded as the superposition of two semi-infinite mining and calculated with Equation 4.
are parameters of the PIM named inflection point offset to the left and inflection point offset to the right.
From the analysis above, the GM to predict the sinking value of the point with abscissa x on the strike main section of the subsidence basin caused by a limited mining coal seam with the length of D 3 is established as Equation 5.
Where L represents the width for calculation of the coal seam and can be calculated by L D S S sin s in are parameters of the GM named mining influence propagation angle, inflection point offsets in the downhill direction and inflection point offsets in the uphill direction.

Parameters analysis
From the analysis above, it can be seen that there needs seven parameters to predict mining subsidence with the GM: . Take a working face in Huainan mining area which has the length of 770m, width of 205m, depth of 556m, thickness of 3m and dip angle of 5° as an example to compare the model curves of the GM and the PIM. The PIM parameters of this working face are as follows: The strike main section predicts the result of using the PIM in combination with the above parameters is shown in curve 1 in Figure 2. And the result of using the GM in combination with the above PIM parameters is shown in curve 2 in Figure 2. By comparing curve 1 and curve 2 in Figure 2, it can be seen that the shapes and trends of the GM and PIM model curves are similar under the same set of parameters, but the calculated specific subsidence values are different. Then adjust the GM parameters to make r r GM PIM  / . 2 5 and again use the GM to predict the subsidence of the strike main section to get curve 3.

Improved cuckoo search algorithm
In order to further study the characteristics of the GM and establish the relationship between the parameters of the GM and the PIM, a suitable algorithm to solve the model parameters is needed. The cuckoo search algorithm (CSA) is a new kind of heuristic bionic group intelligent optimization algorithm proposed by Professor Xinshe Yang and Susan Deb from the University of Cambridge in 2009. The biological behaviors used include parasitic reproduction and Lévy flights mechanism. The CSA has the advantages of clear theory, few parameters, easy to expand, strong global search ability, easy to implement, and so on. It has been deeply studied and widely used by scholars in China and abroad.

Parasitic reproduction behavior
Parasitic reproduction behavior means cuckoos do not build their own nests, they lay their eggs in communal nests or other birds' nests. So mother cuckoos must face the risk that their eggs may be removed from the nests by other birds or the host birds. In the long-term evolution process, some cuckoos have learned some ways to improve the survival rate of their eggs such as pushing the host birds' eggs out of the nests.

Lévy flights
Studies have shown that the foraging routes of many flying animals, including cuckoos, follow the Lévy flights mechanism. The advantage of Lévy flights mechanism in solving optimization problems is that it is combined with high frequency short step movements and occasionally long step movements. The high frequency short step movements make it possible to outcrop the local optimal solution and the occasionally long step movements ensure that the results will not fall into the local optimal solution.

Basic model
According to the above biological background, make the following scenario hypothesis:

•
The total number of the nests is determined; • The cuckoos lay only one egg in each nest; • The one with the best fitness in the nest group will be retained directly to the next generation When the CSA is used in solving optimization problems, the location of each nest represents a set of solutions. After several rounds of elimination, evolution and selection, the optimal individual can be selected out. The elimination strategy of the CSA was random selection under fixed elimination rate. The selection strategy of the CSA was greedy selection, that is, the individuals with better fitness will be retained after the new solutions were obtained and compared with current solutions. The evolution strategy of the CSA adopts two modes: random walking mode and preferred walking mode.
1. Random walking mode. Use Equation 7 when current nests' positions need to be updated Where x i t and x i t1 represent the location of the i-th nest of generation t and Generation t+1 respectively; α represents the step size controller;  is point-to-point multiplication; n is the amount of the nests, that is, the amount of the feasible solutions; L λ ( ) represents the random search path, 3 , where  represents the random step size obtained by Lévy flights.
2. Preferred random walking mode. Mimicking the behavior of a host bird perceives the cuckoo egg and throws it away, then a new nest is needed. The next updating operation is as Equation 8.
Where x j t and x k t represent two random solutions of generation t; r is a random number satisfying the uniform distribution principle between [0,1].

The ICSA
Practice shows that the CSA is not strong in local search ability. This is because each solution is independent of each other in the process of searching for the optimal solution, there is no experience learning from better solutions and no information exchange. Therefore, the roulette selection strategy was adopted instead of the random elimination strategy in the CSA and the average operator and crossover operator was added in the updating process to exchange the information between nests. The specific operations are as follows: The selection strategy.
To select a better individual from the contemporary nests, it is first necessary to establish an individual fitness calculation function. According to the principle of least squares, for the purpose of making the fitting effect of the inverted parameters better, the fitness of solutions are evaluated based on the sum of the squares of the difference between the sinking value of predicted and observed. The fitness function is constructed as Equation 9. The nest with high fitness value is of good quality. Declare a temporary nest X k � for communicating information between current nests. For any individual X t and X t1 , if X X t t = +1 , calculate the step length by α =1 2 / t to generate a new nest X k . If X X t t ≠ +1 , then deal them with three different operators: average operator (take the average of the parameters in the two nests); crossover operator (randomly exchange the parameters in the two nests); and mutation operator (update the nests with Eqs. (8)). The nest X k which has the best fitness can be selected out. Then randomly select a nest from all nests and compare it with X k to keep the one with the better fitness.
Set the total number of the nests as n=200 and each nest represents a set of parameters, which is [ q, tan , , S 1 , S 2 , S 3 , S 4 ] T in the GM. The specific steps of calculating mining subsidence prediction parameters with the ICSA are shown in Figure 3.

Test by simulated working face
A simulated working face is designed so as to test the stability and reliability of the ICSA and compare its accuracy with the CSA and the GA,.
Its length, width, depth, dip angle and thickness are designed as 900m, 300m, 400m, 0° and 3m. Its GM parameters are designed as is shown in Table 1. There set up an observation line above both the strike main section and inclined main section as is shown in Figure 4. The spacing and total number of the observation points are 30m and 88. Take the value calculated by the GM with designed parameters as true subsidence value. Then add random errors with median error of 5m, 10mm and 15mm into the true subsidence value to generate three sets of designed observations. In order to test the stability of the ICSA, calculate 100 sets of the GM parameters of the designed working face with the above designed observations and the comparison of the parameter inversion effect of the ICSA while the designed observation data has different degrees of median error are shown in Table 2. In order to test the precision of the ICSA, use the GA, the CSA and the ICSA to calculate 50 sets of the GM parameters of the designed working face with the above designed observations. Take "q" as an example, the error comparison between true value and values inversed by the GA, CSA and ICSA are shown in Table 3.    We can draw the following conclusions from Table 2 and Table 3: 1. the relative errors of q GM , tan GM  and  GM inversed by the ICSA are all less than 6%, and the errors of S are all less than 10mm. Indicating that the errors of the parameters inversed by the ICSA are all within the normal range. The average values of the errors of q GM , tan GM  and  GM are all less than 1, indicating that the results calculated by the ICSA are stable and reliable.
2. as the medium error of the random error in the observation data increases, the medium error of the inversed parameters has a tendency to increase. Indicating that the quality of the observation value will directly affect the accuracy of the parameters inversed by the ICSA.
3. the average error of q GM inversed by the ICSA is 0.035, and its relative error is 5.0%, which is similar to the GA and significantly higher than the CSA, indicating that the stability of the ICSA is higher than the CSA and  comparable to the GA. Compared with the other two methods, the error in the fitting of the parameters calculated by the ICSA is the smallest, indicating that the reliability of the parameters calculated by the ICSA is the highest.

Case analysis
(1) The working face 1111 (1) of the Zhujidong coal mine in Huainan Panji mining area has a length of 1584m, width of 220m, average depth of 920m, average thickness of 1.8m, and average dip angle of 3°. It is buried under the loose layer with an average thickness of 278m. Its coal mining method and roof management method are the comprehensive mechanized method with a full height mining at one time and the total collapse method respectively.
Using the ICSA combined with the last observation data of the working face to calculate the PIM parameters and the GM parameters. The calculated model parameters are shown in Table 10. The comparison between the observations and the subsidence value fitted by the GM and the PIM are shown in Figure 5 and Table 4.  (2) The working face 12326 of the Gubei coal mine in Huainan Xieqiao mining area has the length of 770m, width of 205m, average depth of 537m, average thickness of 2.55m, and average dip angle of 5°. It is buried under the loose layer with an average thickness of 448m. Its coal mining method and roof management method are the comprehensive mechanized method with a full height mining at one time and the total collapse method respectively.
Using the ICSA combined with the last observation data of the working face to calculate the PIM parameters and the GM parameters. The calculated model parameters are shown in Table 10. The comparison between the observations and the subsidence value fitted by the GM and the PIM is shown in Figure 6 and Table 5.  (3) The working face 1613 (1) of the Guqiao coal Mine in Huainan Xieqiao mining area has a length of 1528m, width of 251m, average depth of 668m, average thickness of 2.8m, and average dip angle of 3°. It is buried under the loose layer with an average thickness of 420m. Its coal mining method and roof management method are the comprehensive mechanized method with a full height mining at one time and the total collapse method respectively.
Using the ICSA combined with the last observation data of the working face to calculate the PIM parameters and the GM parameters. The calculated model parameters are shown in Table 10. The comparison between the observations and the subsidence value fitted by the GM and the PIM is shown in Figure 7 and Table 6.  (4) The working face 1414(1) of the Guqiao coal Mine in Huainan Xieqiao mining area has a length of 2120m, width of 251m, average depth of 735m, average thickness of 2.7m, and average dip angle of 5°. It is buried under the loose layer with an average thickness of 411m. Its coal mining method and roof management method are the comprehensive mechanized method with a full height mining at one time and the total collapse method respectively.
Using the ICSA combined with the last observation data of the working face to calculate the PIM parameters and the GM parameters. The calculated model parameters are shown in Table 10. The comparison between the observations and the subsidence value fitted by the GM and the PIM is shown in Figure 8 and Table 7.  (5) The working face 1252(1) of the Panyidong coal mine in Huainan Panji mining area has the length of 1148m, width of 262m, average depth of 802m, average thickness of 2.7m, and average dip angle of 6°. It is buried under the loose layer with an average thickness of 165m. Its coal mining method and roof management method are the comprehensive mechanized method with a full height mining at one time and the total collapse method respectively.
Using the ICSA combined with the last observation data of the working face to calculate the PIM parameters and the GM parameters. The calculated model parameters are shown in Table 10. The comparison between the observations and the subsidence value fitted by the GM and the PIM is shown in Figure 9 and Table 8.  9. Comparison of the fitting effect of the GM and PIM for 1252(1) working face (6) The working face 1312 (1) of Gubei coal mine in Huainan Xieqiao mining area has the length of 629m, width of 205m, average depth of 528m, average thickness of 3.6m, and average dip angle of 5°. It is buried under the loose layer with an average thickness of 439.7m. Its coal mining method and roof management method are the comprehensive mechanized method with a full height mining at one time and the total collapse method respectively.
Using the ICSA combined with the last observation data of the working face to calculate the PIM parameters and the GM parameters. The calculated model parameters are shown in Table 10. The comparison between the observations and the subsidence value fitted by the GM and the PIM is shown in Figure 10 and Table 9.   The six working faces selected in this article all buried under a thick loose layer. We can draw the conclusion from the above calculation results that the GM fitted better than the PIM at the boundary of the subsidence basin. Taking 1613 (1) as an example, the median error at the boundary of the subsidence basin fitted by the PIM is 63.7, while fitted by the GM is 29.4, which is a reduction of 34.3; The relative error at the boundary of the subsidence basin fitted by the PIM is 57.1%, while fitted by the GM is 23.2 %, which is a reduction of 33.9%. The fitting effect of the GM at the boundary of the subsidence basin is significantly optimized compared with the PIM in predicting the subsidence value caused by mining under thick loose layers. By comparing the GM parameters and the PIM parameters in Table 9, it is confirmed that the values of the two models parameters are basically the same except for � r GM and r PIM . By comparing their values in Table 9, it can be seen that r PIM is about 2.5 times that of � r GM .

Results and discussion
(1) Through comparison, it is proved that the value at the boundary of the GM's function curve is larger than that of the PIM. Which makes the GM more suitable than the PIM for predicting surface subsidence caused by mining under thick loose layers.
(2) In this paper, the ICSA was established by adding operations to increase the information exchange between individuals on the basis of the CSA. And through a simulated working face, it is proved that the stability and accuracy of the ICSA can meet the requirements of coal mining subsidence prediction parameters inversion.
(3) 6 working faces buried under thick loose layers in Huainan mining area were selected in this paper to compare the fitting effect of the PIM and the GM using the ICSA. The results show that the GM has a better fitting result than the PIM at the boundary of the subsidence basin. Further confirmed that the GM is more suitable than the PIM for predicting surface subsidence caused by mining under thick loose layers. And the parameters of the GM and the PIM can take the same value except for r r GM PIM  / . 2 5. This result helps the GM to be used in actual engineering