Numerical study of reasonable cycle step length for longwall top-coal drawing in extra-thick coal seams based on the particle-block element coupling approach

: The cycle step length (CSL) is a significant parameter for longwall top-coal drawing technology that remarkably affects the top-coal recovery rate and the rock-mixing rate, especially for extra-thick coal seams. In this study, a particle-block element coupling approach is performed to investigate a reasonable CSL for extra-thick coal seams. By comparing this approach to the Bergmark-Roos analytical result, the proposed numerical model is verified, showing good performance in modeling top-coal caving. A 2-D numerical model of hydraulic support considering the mechanical behavior of the legs is established, which can be used for modeling the interaction between the hydraulic support and the top coals during the top-coal drawing process. The top-coal recovery rate, the top-coal drawing body shape, and the evolution characteristics of the coal-rock interface under different CSL conditions are compared. In addition, the mechanism of the lost top coal affected by the CSL is revealed. The results show that the CSL of top-coal drawing has a significant effect on the morphology of the coal-rock interface and the mutual invasion of coal and rock, which is the primary reason for coal loss and further affects the top-coal recovery rate and the rock-mixing rate. It is suggested that the CSL should be 0.8 m when the top-coal thickness is 12 m.

. Numerical simulation has also developed rapidly into an effective method for studying LTCC. Khanal et al. (2011Khanal et al. ( , 2014 established a numerical LTCC model with COSFLOW software and analyzed the support collapse, top coal failure mechanism, roof caving mechanism, support stress, and vertical stress on the caving effect and finally discussed the feasibility of using LTCC technology in the Adriyala coal mine. Simsir F et al. (2008) simulated the coal production of an LTCC working face using ARENA software and found that the one-per-two reamer method has a higher top-coal recovery rate. In addition, the finite difference method was used to analyze the mechanism of top-coal caving in the Omerler coal mine (Yasitli and Because of the complex environment behind hydraulic supports, it is difficult to observe the top-coal caving process above the hydraulic support. A numerical simulation is thus a fast and effective way to simulate the top-coal caving process. The particle element method is the major simulation method used for studying top-coal drawing, but modeling the interaction between the hydraulic support and the top coal during the drawing process presents a challenge. In this study, the particle-block element coupling approach is used to simulate the longwall top-coal drawing process. The block mesh model and the constitutive law of hydraulic support are proposed, and the interaction principle between particles and block elements is introduced. The process of top-coal drawing with a single coal drawing opening is modeled and verified with the Bergmark-Roos analytical model. Then, the correlation between the top-coal recovery rate and the shape of the coal-rock interface is analyzed, and based on this result, a reasonable CSL for extra-thick coal seams is determined.

The interaction principle between particle and block element
In a numerical simulation method of top-coal mining, the flow of particles from aggregates broken into loose and falling particles can well simulate the top-coal drawing process and the coal-rock interface state. In addition, a hydraulic support numerical model is added to the simulation process. The block element can accurately simulate and analyze the stress state of the coal rock and hydraulic support so that the simulation result of the top-coal drawing process is closer to the actual situation of the project. The calculation principle of particle and block coupling used in this paper is as The internal mechanical calculation of the block element adopts the explicit solution method to calculate the finite element, including the calculation of the element node force and the calculation of the element node motion.
The calculation formula of the resultant force on the element node is as follows: where F is the resultant force of the element node, Fe is the external force of the element node, FD is the deformation force of the element node, and FC is the damping force of the element node.
The calculation formula of node movement is as follows: where is the acceleration of the node of the calculation element; is the velocity of the node of the calculation element; ∆ is the displacement increment of the node of the calculation element; is the full amount of the node displacement of the calculation element; is the mass of the node of the calculation element; and ∆ is the node of the calculation element time step.
Using the incremental calculation method to calculate the stress of the element and the node deformation of the element gives the following: where The linked bar model can be used to simulate the continuous medium properties between the particles. A connecting rod model is established between the two particle elements, as shown in Figure 1. It is used to simulate the process from the initial intact state to crushing of the top-coal drawing. The link bar model is regarded as a rectangle, through which the contact force or cohesive force between the two particle elements can be calculated. The long side of the rectangle is the sum of the radii of the two particles, while the short side is equal to the diameter of the smaller particles. Based on the connecting rod model, the relationship between the two particle elements is the surface contact, and the equivalent contact area AC is the projected area of the smaller particles. The formula for calculating the contact force between the discrete elements of the particles is as follows: where ∆ and ∆ are the incremental difference of the normal displacement and the incremental difference of the tangential displacement between two discrete elements of particles in contact with each other.
The formula for calculating the contact torque between discrete elements of particles is as follows: where and are the torque and bending moment between the discrete elements of the particles, I and J are the moment of inertia and the moment of inertia between the contact surfaces of the discrete elements of the particles, and ∆ and ∆ are the incremental differences between the torsion and bending angles between the discrete elements of the particles.
where 1 and 2 are the radius values of the two particles contacting each other between the discrete elements of the particles, and is the contact area between the discrete elements of the particles. The contact stiffness between the discrete elements of the particles can be derived from the elastic modulus and shear modulus of the particles in contact with each other: where and are the normal and tangential stiffness ̅ between the discrete elements of the particles in contact with each other, and ̅ is the average elastic modulus and the shear modulus of the two discrete elements of the particles in contact with each other.
According to the Mohr-Coulomb criterion and the maximum tensile stress criterion, the contact force calculation formula is as follows: The contact judgment condition between the two particle elements is formula (1.9).
If any one of the inequalities in the formula is satisfied, the contact between the particles will no longer transmit torque.
where = ( 1 + 2 )/2， ， ，and are the tensile strength, cohesion and internal friction angle, respectively, and I is the moment of inertia.
The calculation of the torque on the particle discrete element is as follows: where 1 and 2 are the vectors of the rotational angular velocity of particle discrete elements 1 and 2; 1 and 2 are the vectors of the relative positions of particle discrete elements 1 and 2 to the contact point (from the particle centroid to the contact point); and 1 and 2 are the translational velocity vectors of the centroid of particle discrete elements 1 and 2.
where 1 and 2 are the torque on particle discrete elements 1 and 2, respectively, and ( ) is the contact force of the particle discrete element in the global coordinate system.
The core of the coupling calculation of the block element and the particle discrete element is the logical judgment of the mutual contact between the block element and the particle discrete element. In the two-dimensional numerical calculation, the method for judging the contact between the block element and the particle discrete element is the body center of the particle element. To judge the relative position of contact with the edge of the block element, the contact between the body center of the particle element and the edge of the block element must also satisfy that the distance from the body center of the particle discrete element to the boundary edge of the block element is less than or equal to the radius of the particle discrete element. That is ≤ , and the projection point of the body center of the discrete element of the particle on the boundary edge of the block element is inside the edge of the block element, that is, ≤ , ≤ can be established. Once the particle discrete element and the edge of a boundary block element have established a contact relationship, the normal spring and tangential spring that contact each other between the block element and the particle discrete element are automatically created, and the block element contacts the particle discrete element. The interpolation coefficient of point k will be automatically calculated by the following formula: where is the relative position vector of the particle from the element center k and the block element edge i, i and j are the two endpoints of the block edge, djk is the distance between point j and point k, dij is the difference between point i and point j, dik is the distance between point i and point k, and n is the normal vector outside at the edge of the block.
According to the alternating cycle calculation process of the theoretical formula (1~12), the explicit solution process for the finite element, discrete element, particle, and block coupling can be realized.

The constitutive law of hydraulic support
In previous simulation methods for top-coal drawing, most scholars used PFC to simulate the top-coal drawing process. To more realistically simulate the top-coal drawing process at the top-coal drawing working face, the constitutive model of the top-coal hydraulic support is put into the CDEM simulation. According to the coupling principle of block and particle introduced in Section 2.1, the interaction between the coal gangue particles and the hydraulic support for the top-coal drawing and the change process of the working resistance of the hydraulic support during the top-coal drawing process are simulated to more truly simulate the on-site coal caving process. The test points are arranged between the inner and outer columns of the hydraulic support, as shown in Figure 2(a). The red point is the measuring point of the working resistance of the hydraulic support for top-coal drawing. The hydraulic support model adopts the constitutive model, which can represent the relationship between the column shrinkage and the working resistance. The calculation formula is as follows: where P is the working resistance of the support; P0 is the setting force of the support; K is the hydraulic stiffness; ΔS is the shrinkage of the column.
The numerical constitutive curve of the hydraulic support is shown in Figure 2 Figure 2(b) shows that the constitutive structure of the stent in the simulation has two stages. In the first stage, the support resistance from the initial support force P0 to the working resistance P1 with stiffness K linearly increases as the column shrinks down. The safety valve opens in the second stage, and the working resistance of the support will not change as the column continues to shrink. As the working face advances, the hydraulic support is recycled from the above process.
The incremental displacement (m) The working resistance (MPa) The hydraulic stiffness , K Softening , K 1 , K 2 , K 3 Rated working resistance , P 1 Initial support force, P 0 According to Figure 3(a), the distance from the origin of the coordinates to any point on the opening on the drawing is as follows: where Putting equations (15) and (16) into (14), we obtain: where r1 is the distance from the center of the drawing opening to the origin of the polar coordinates; r2 is the distance from any point on the drawing opening to the origin of the polar coordinates; is the angular coordinate of the particle; is the maximum allowable displacement boundary angle. The force is equal to the angle of the particle weight; when = 0, rmax is the maximum distance from the origin of polar coordinates to the boundary of the drawing body; and D is the size of the opening in the figure.

Comparison
According According to the ellipsoid theoretical simulation equation of formulas (16) and (17), the single coal drawing port model shows that the thickness of the top coal is 1 =12 m, and the size of the coal drawing opening is d = 1.75 m. From the numerical simulation of a single coal drawing opening, =30° can be obtained, and equation (17) is simplified to equation (18): According to formula (18), the theoretical model of the shape of the top-coal drawing body parallel to the cross-section of the working face is shown in Figure 3 Table 2 and Figure 4.

Numerical simulation scheme
To study the influence of different CSLs on top-coal recovery and the evolution of

Determination of the best coal caving step
The drawn amount of top-coal drawing in each drawing cycle under different CSL conditions is shown in Figure 8 below. It can be seen that, regardless of the adopted coal drawing CSL, the first drawn amount is the largest, and the drawn amount increases with the increase of coal drawing CSL. Figure 9 Advance distance of working face (m)  Table 3.
The top-coal recovery rates from 11 different drawing CSLs are compared and shown in Figure 10 below.

Conclusion
In this study, a particle-block element coupling approach is implemented to investigate a reasonable cycle step length of top-coal drawing for extra-thick coal seams. The block mesh model and the constitutive law of hydraulic support are proposed, and the interaction principle between particles and block elements is introduced. The correlation between the top-coal recovery rate and the shape of the coal-rock interface is analyzed, and based on this result, a reasonable CSL for extra-thick coal seams is determined. The main conclusions are as follows: (1) By comparing the numerical simulation result to the Bergmark-Roos analytical result, the proposed numerical model is verified and shows good performance in modeling top-coal caving. A 2-D numerical model of hydraulic support considering the mechanical behavior of the legs is established, which can be used for modeling the interactions between hydraulic support and top coals during the top-coal drawing process.
(2) The distribution of the lost coal in gob and its mechanism are clarified. During the advance process of top-coal drawing, with the rock reaching the coal drawing opening and with periodic floating, the top-coal drawing body is affected by the tail beam of the hydraulic support, which takes the shape of a cutting variation ellipsoid. In addition, the displacement angle of the ellipsoid to one side of the gob is different with different CSLs. If the CSL is too long, the upper rock will arrive at the coal drawing opening before the coal behind the hydraulic support. If the coal drawing CSL is too short, the rock in the gob will rush into the coal drawing opening in advance. As a result, the loss of top coal is periodically inclined to one side of the gob. In each drawing cycle, the top coal flows regularly twice. The first flow is a large range of top-coal falling after the support movement, and the other flow is the process of top coal upper rock flowing with the top coal at the coal drawing opening.
(3) Combined with the actual geological situation of the no. 8222 working face, the top-coal recovery rate, the shape of the top-coal drawing body, and the evolution characteristics of the coal-rock interface under different CSL conditions are compared. In addition, the mechanism of the lost top coal affected by the CSL is revealed. The results show that the CSL of top-coal drawing has a significant effect on the morphology of the coal-rock interface and the mutual invasion of coal and rock, which is also the primary reason for coal loss and further affects the top-coal recovery rate and the rock-mixing rate. It is suggested that the CSL should be 0.8 m when the top-coal thickness is 12 m.

Declarations
We declare that we do not have any commercial or associative interest representing a conflict of interest in connection with the paper submitted.