Analysis of the Influence of Fault Fracture Zone on Mining Response Based on FDM-DEM Coupling

Fault slip will cause a change in mining stress at the longwall face, which will cause adverse e ﬀ ects. In this study, on the basis of Fast Lagrangian Analysis of Continua in 3 Dimensions (FLAC3D) and Particle Flow Code in 3 Dimensions (PFC3D), the sliding of the fault fracture zone and its impact on the longwall working face were analyzed. The rock mass of the fault fracture zone with a certain thickness was constructed using rigid random model particles. The coupling between the wall element of PFC3D and the zone element of the continuous medium in the software was used to realize the transmission of force and displacement, and the interaction between the fault fracture zone and the working face was studied. The in ﬂ uence of the slip of di ﬀ erent fault zone positions on the fault and working face was also explored using the method of externally disturbing the fracture zone of the fault. The numerical results showed that as the distance between the fault and the working face continued to decrease, the peak stress concentration in front of the longwall face ﬁ rst increased, then decreased, and gradually shifted to the vicinity of the fault zone. The stress mutation and fault slip occurred within a certain distance of the longwall face from the fault. When fault slip activation begins, the stress near the fault zone showed a sudden change of varying degrees with the advancement of the longwall face. This sudden change was caused by the in ﬂ uence of mining activities from the activation distance of the rock in the fault fracture zone and the rolling extrusion of the rock mass in the fracture zone. When the fault zone closer to the working face was disturbed, the in ﬂ uence on the fault zone and longwall face was greater. When the fault zone near the coal seam was disturbed, the rock mass and working face near the fault zone brought di ﬀ erent degrees of dynamic responses, which were mostly instantaneous and had high frequency and amplitude. The results of this research could help in the mining of longwall mining a ﬀ ected by fault zones and have a certain guiding role in coal mining before crossing faults.


Introduction
Fault is a type of geological structure often encountered in longwall mining. Underground coal mining will cause redistribution of surrounding in situ stress [1]. Changes in in situ stress may activate nearby faults, resulting in the slip and instability of the faults, and the activity of the faults will affect and bring certain difficulties to the mining of coal [2,3]. From past experience, the fault slip often leads to coal mine disasters, such as rock bursts [4], water penetration [5], and coal and gas outbursts [6][7][8][9]. Therefore, the study of fault slip in coal mining can help in thoroughly under-standing the characteristics of fault slip in underground mines to prevent the risks caused by fault slip.
In previous studies, faults were often treated as plane of weakness [10][11][12]. The mining subsidence fault layer blocks the stress and energy transfer caused by mining, but in actual geological conditions, faults often appear in the form of fault zones with a certain thickness [13][14][15][16]. Wei et al. [17] used the Fast Lagrangian Analysis of Continua in 3 Dimensions (FLAC3D) model to analyze the source mechanism of the fault slip behavior. The postfailure process of the fault slip was studied by constructing a zero-thickness interface element to simulate the fault, and they concluded that the fault slip shows distinct trends in different positions of the fault. To investigate the sudden failure of faults, Wang et al. [18] performed a physical similarity simulation of a double-fault structure. The simulation results suggested that changes in strain rates from relatively low to abrupt changes may be identified as precursors to fault slip. To study the deformation characteristics of the surrounding rock in the fault zone of a coal mine, Kang et al. [19] established a 3D model of the fault roadway and used special support to treat the roadway and grouting for the fracture, which provided favorable conditions for roadway stability.
Fault planes are often presented in the form of interface element in numerical calculations [20], and fault zones are often replaced with weak interlayer zones; however, most of actual fault zones are discontinuous fragments [21]. Numerical methods are widely used in longwall mining research, and the coupling method of finite difference method (FDM) and discrete element method (DEM) has been well used in geotechnical engineering [22][23][24][25][26]. Li et al. [27] established a 3D coupling model considering both DEM and FDM, simulated ballast particles with DEM considering rolling resistance model, realized the coupling process between DEM and FDM with interface element, and verified the numerical model under impact load through field test. Feng et al. [28] used Particle Flow Code in 3 Dimensions (PFC3D) and FLAC3D to simulate the landslide process. The seismic signal characteristics of the Xiaolin landslide in 2009 were studied and compared with the broadband seismic station data.
The fault slip is caused by the mining activities or a change in tectonic stress. The specific breaking position of the fault is not clear. In this study, the coupling method of FDM and DEM was used to establish a 3D geological fault model. The mining response under the influence of the fault fracture zone was simulated, and the method of external disturbance fault slip was used to deal with the fault fracture zone, to study the response under quasi-static and dynamic conditions, and to provide certain guidance for coal mining under fault conditions.

Study Area
Jisan Coal Mine is located in southwest of Jining City, Shandong Province, China. As shown in Figure 1(a), the mine covers an area of 110 km 2 . It is estimated that the recoverable coal reserves are 526 million tons. The initial design of the mine type is 5 million tons per year. The output once reached 8 million to 10 million tons per year. The 6303 working face is arranged in the west of the 6300, 6301, and 6302 working faces in the sixth mining area (Figure 1(c)). The designed inclined length is 239.8 m, and the strike length is 2057.8 m. The working face return airway is a gob-side tunnel with a total length of 2114 m and a working face buried depth of 669~711 m. The 6303 return airway is adjacent to the 6302 gob, and the 6302 gob is connected with the 6303 working face return airway. Narrow coal pillars of 2~3 m are set between them, the average thickness of the coal seam is 4.75 m, and the coal seam is a near-horizontal coal seam. The coal seam column chart is shown in Figure 1(b). During the advancing process of the working face, rock burst events occurred many times. During the excavation of the return airway along the trough, a total of three faults SF62, SF28, and SF3 were passed, and the fault drops were 1 m, 1.6 m, and 0.5 m, respectively. zone was also replaced with a continuum weakened zone with a vertical height of 14 m above the lower boundary of the model, and the middle part was modeled by broken rock mass. The strata and weakened zone in the model were modeled with finite difference continuous elements. The physical and mechanical parameters of coal and rock were selected with reference to the C8-9 drilling data of the Jisan Coal Mine [32]. The specific parameters are presented in Table 1, and the 3D size of the numerical model is shown in Figure 2.

Modeling and Parameter Selection for the Fracture Zone.
Irregularly shaped gravel was used in the simulation of the fault fracture zone proposed in the numerical model to simulate the discrete properties of the material. In consideration of the force characteristics of the fault zone, the gravel in the crushed zone was assumed to be rigid; that is, it could not be broken again, and a single gravel was regarded as a rigid block that would not fail. PFC3D introduces the method of nonspherical particles. Therefore, the rigid block in the PFC3D software can be used to simulate the fault fracture zone [33]. A single modeling of the fault fracture zone was based on the 3D particle template drawn by CAD, and several random templates were generated using the random film variation coefficient method. The template was imported into the PFC3D software to generate first rigid block particles and then the fault fracture zone. The particle mechanical parameters of the fault fracture zone in the model were selected in accordance with those of gravel in the literature [34]. The computational efficiency for rigid blocks was maximized because of the contact detection and resolution of the Gilbert-Johnson-Keerthi algorithm [35], which assumes only one contact between two protruding objects. The simulation of the fault fracture zone is shown in Figure 3. The contact behavior within the fracture zone of the fault was represented using a linear contact bond model that provides a linear spring acting in parallel with two fractured rocks [36], and the specific micromechanical parameters are shown in Table 2. The coupling calculation scheme adopts the wall-element cou-pling method in FLAC3D [37]. The principle of coupling is to use the contact force and moment of the wall to determine the equivalent force system at the contact point. Force is passed to the grid node along with stiffness, and force and stiffness are logically checked and updated. The updated force and stiffness have changed the geometric parameters and cell stiffness of the solid element to ensure the stability of the calculation process. With reference to the description of the FLAC3D software, in accordance with the parameters of the near elements, the  3 Geofluids parameters of the interface could be based on the equivalent stiffness of the adjacent zone elements, that is, The selection was conducted. K is the bulk modulus, G is the shear modulus, and Δz min is the minimum with of adjacent regions in the normal direction. The k n and k s of interface parameters is 0.66 MPa, and the value of internal friction angle is 15°.

Calculation and Mining Influence Analysis under Static Load Conditions
The coal mining method for the working face was the longwall mining method, and the calculation under the static load condition of the model was divided into two parts. First, in situ stress balance of the model was performed. After the in situ stress balance was completed, the overall displacement of the model was initialized to zero. Then, excavation calculation of the coal seam was carried out. The model was excavated from the right boundary, and the farthest distance from the fault plane was 200 m. To observe the stress and fault changes when the working face approaches the fault, the excavation of the model was divided into two stages. In the first stage, the model excavation was 20 m per excavation step within the range of 200-60 m longwall face from the fault. In the second stage, the model excavation was 10 m per excavation step within the range of 60-10 m longwall face from the fault. In the process of coal seam mining, the working face was constantly approaching the fault. The specific excavation design is shown in Figure 4.

Coal Seam and Working Face Stress of the Longwall
Face. The vertical stress change and nephogram in the coal seam during the mining process in the first stage are shown in Figure 5. During excavation stage 1, due to the distance from the fault zone, the coal seam and working face were less affected by the fault. During the mining process of the longwall face, the vertical stress in the coal seam first increased and then decreased from the longwall face to the fault. The stress peak appeared within a certain distance in front of the longwall face. With the continuous advancement of the longwall face, the stress peak in front of the working face continued to rise. When the coal seam was mined for 20 m, the stress peak of the coal seam was 29.3 MPa at 7 m in front of the coal wall of the working face. When the coal seam was mined 80 m away from the fault plane, the maximum stress peak at 10 m in front of the working face was 63.4 MPa. As the distance between the working face and the fault decreased gradually, the change trend of stress was small in the range from the working face to the peak stress. In the range from the peak to the fault, the downward trend of stress gradually increased, which was caused by the gradual increase in the stress peak in front of the working face with the continuous advancement of the working face.
The stress change and nephogram in the coal seam during the second stage of coal seam mining are shown in Figure 6, which depicts that the stress change in front of the working face first increased and then decreased, and the stress peak was in front of the working face. However, the peak stress in front of the working face gradually decreased as the distance between the working face and the fault plane decreased. When the longwall face was 60 m away from the fault, the stress peak at 8.5 m in front of the working face was 57.6 MPa; when the working face was 10 m away from the fault, the stress peak at 6.5 m in front of the working face was 21.6 MPa. The stress nephogram illustrates that with the gradual reduction in the distance between the working face and the fault, the stress peak in front of the working face gradually shifted to the vicinity of the fault zone, mainly because the fault zone was a fracture zone. The mining of the working face caused a stress change in the fault fracture zone, slip of the fault zone, and stress release of the fracture zone. Figure 7 shows the stress peak and the change in its position from the working face during the advancing process of the whole longwall face. When the longwall face was within 200-80 m from the fault, the peak stress gradually increased; when it was within 80-20 m from the fault, the peak stress gradually decreased. As the working face approached the fault, the stress was affected by the fault and tended to rise again. The position of the peak stress in front of the longwall working face gradually increased within 200-100 m, then decreased within 100-30 m, and increased again within 30-10 m of the working face from the fault. That is, when the working face advanced to 60 m or more away from the fault, the stress peak and its position from the working face were affected by the conventional mine pressure. When the working face advanced to 60 m away from the fault and gradually approached the fault, the stress peak of the working face was gradually affected by the fault slip.

Displacement
Response Analysis. Figure 9 shows the vertical displacement of the eight monitoring points and the displacement difference between adjacent points on both sides of the fadult zone. The fault started from the activation distance, and the displacement of the monitoring points on both sides began to change.

Static Analysis of the External Disturbance Fault Fracture Zone
External disturbance is the main cause of fault slip [34]. Under the condition of fault, underground coal mining may cause fault slip, but the specific slip position of the fault cannot be determined. In previous studies, experiments or numerical simulations were often used to perturb the fault [38,39]. Observing the slip characteristics of the fault, this study conducted external disturbances in the fault fracture zone. The specific disturbance positions are shown in Figure 10.  Figure 11 is the vertical stress variation diagram of four rock stratum monitoring points in the footwall of the fault under external disturbance. Figure 11(a) is the stress variation of point A, which was far away from the coal seam vertically. Owing to the distance from the coal seam and being inconsiderably affected by the mining of the longwall working face, the change in vertical stress under external disturbance was insignificant. Figure 11(b) shows the vertical stress change of monitoring point B of the rock stratum at a vertical distance of 30 m above the coal seam. This layer of rock was closer to DP4, so DP4 had a more significant influence on the vertical stress of the monitoring point. When the longwall working face was mined to 40 m from the fault, the stress decreased from 16.9 MPa to 14.4 MPa; when the distance was 30 m, the stress decreased from 15.7 MPa to 11.7 MPa. Figure 11(c) illustrates the vertical stress change at monitoring point C (5 m above the coal seam). DP3 had a significant impact on roof rock layer C. When the longwall working face was 60 m away from the fault, the stress during the disturbance dropped by 5 MPa. When the working face was 10 m away from the fault, the stress dropped from 40.3 MPa, which was not disturbed, to 28.2 MPa. The rest of the disturbance points had no significant effect on the vertical stress of the monitoring point. 9 Geofluids generated at DP1 was 8.2 MPa, and the maximum drop of the vertical stress at DP2 was 9.3 MPa compared with the undisturbed one. DP1, DP2, and DP3 exerted a significant impact on the rock formation disturbance near the coal seam. When the slip position of the fault was closer to the coal seam, it might have a greater impact on the stress change near the mining area. Figure 12 is the vertical stress nephogram of the coal seam and the vicinity of the fault under three external disturbances. The previous calculation indicated that the closer the fault disturbance was to the fault, the more significant the impact would be. DP4, which was far from the vertical height of the coal seam, had little impact on the coal seam, so it had no value. The cloud map was drawn 30 m from the longwall face to 10 m from the fault. Figures 12(a), 12(d), 12(i), and 12(l) are the vertical stress nephograms of the undisturbed and disturbed positions DP1, DP2, and DP3 when the local fault of the longwall working face was 30 m, respectively. Different disturbance points had diverse effects on the stress. In particular, DP1 exerted a minimal effect on the stress concentration in front of the longwall working face, but the overall stress was reduced by 0.2 MPa, which significantly changed the stress under the coal seam near the fault. DP2 slowed down the stress concentration in front of the longwall working face, which increased the stress of the rock mass near the fault with the hanging wall of the fault parallel to the coal seam to a certain extent. DP3 also alleviated the stress concentration in front of the coal wall of the working face, and the position of this stress concentration gradually moved up to the rock layer of the coal seam roof instead of existing in the coal seam in front. With the continuous advancement of the working face, the distance between the longwall face and the fault was gradually shortened. Under the influence of the fault zone, the stress concentration in front of the coal seam of the working face was reduced, released, and transferred to the fault, and the stress concentration near the fault zone was gradually increased. When the longwall working face was far away from the fault, the stress of the strata in the hanging wall of the fault was higher than that in the footwall. It was also gradually transferred from the hanging wall of the fault to the foot wall of the fault. The reason was that, with the advancement of the longwall face, the size of the complete rock mass between the coal wall and the fault gradually decreased, resulting in a continuous increase in stress and aggravation of stress concentration.

Dynamic Response and Seismic Signal Characteristics of Externally Disturbed Fault Zones
The slip of a fault may lead to the occurrence of instantaneous vibrations, such as mine shocks and rock bursts [40,41]. Dynamic analysis of the model under different disturbance positions was carried out, and the conditions that may cause stress changes were calculated with reference to the above results. The model with a distance of 10 m from the fault on the longwall working face was perturbed, and DP1 and DP2, which were closer to the coal seam, were selected for calculation. The influence of the model itself and the factors of the disturbance of the broken zone was disregarded, and dynamic boundary conditions were applied to the model after the broken particles were deleted. No dynamic load was applied to the model, and the acceleration near the fault and the working face was monitored within the dynamic monitoring time of 1 s. Wavelet transform is the earliest time-frequency signal processing method based on wavelet proposed by Jean Morlet. Wavelet transform can process the signal in time and frequency, accurately grasp the overall characteristics and local details of the signal, and accurately describe the characteristics of nonstationary signals in time-frequency.
Empirical mode decomposition (EMD) is an adaptive time-frequency decomposition method for signal processing proposed by Huang et al. [42], especially for the processing of nonlinear and nonstationary signals. This method processes the signal according to the time characteristics of the signal itself, without setting the odd function in advance. EMD method has been widely used in different engineering fields. The key of EMD is to decompose the complex signal into finite intrinsic mode functions (IMF), and the decom-posed IMF component is the characteristic signal that characterizes the original signal at different time scales. This method is intuitive, direct, and adaptive. Figure 13 shows the vertical acceleration, wavelet decomposition time-frequency, and empirical mode decomposition (EMD) maps of the monitoring working face and the vicinity of the fault when the longwall working face was 10 m away from the fault under the disturbance of DP1. Figure 13(a) is the vertical acceleration of the monitoring point near the fault plane. According to the data in the figure, when the fault was disturbed, a large acceleration amplitude occurred within 0.2 s. Therefore, the sliding impact of rock mass in the fault fracture zone had a certain impact on the rock mass, which was instantaneous. In the period after 0.2 s, the acceleration amplitude was small and stable. Figures 13(b)  11 Geofluids the duration of time, and the high-frequency part of the vibration was mainly concentrated before 0.5 s. The frequency of the vibration also decreased gradually with time. Given that the working face was affected by the overlying rock layer and the fault zone at the same time, the continu-ous vibration time of the longwall working face was longer than that near the fault zone. Figure 14 shows the vertical acceleration, wavelet decomposition time-frequency, and EMD maps of the monitored working face and the vicinity of the fault when the     13 Geofluids instantaneous vibration when the longwall face approaches the fault, which will affect the mining of coal seams.

Conclusions
A FLAC3D-PFC3D-based coupled geological model of continuous discrete media was established, and the fault zone was modeled with discrete media to analyze the influence of the fracture zone of a fault and the sliding mechanism of the fault on longwall mining. The results showed that in the process of the longwall working face approaching the fault zone, the stress peak in front of the longwall face first increased, then decreased, and finally transferred to the fault fracture zone. Owing to the activation distance of the fault, the monitoring points on both sides of the fault had different degrees of stress mutation. This stress mutation mainly occurred in the rock stratum within 0-30 m of the vertical distance of the coal seam. The rock stratum higher from the vertical distance of the coal seam had greater displacement affected by mining, but no stress mutation occurred. On the contrary, the rock stratum of the coal seam floor was less affected by the fault zone and mining.
Four positions were selected to disturb the fault zone externally to judge the influence of the fault slip at different 14 Geofluids positions on the working face and the rock formations near the fault. The numerical results showed that the fault position slip close to the coal seam would have a significant impact on the working face and the fault. This impact was mainly reflected by the sudden drop in stress in the rock stratum near the fault zone, which might cause sudden seismic events. The dynamic response of the disturbanceaffected fault zone was monitored. The disturbanceaffected faults would produce short-time high-frequency vibrations with large amplitudes.

Data Availability
All data, models, or code generated or used during the study are available from the corresponding author by request (list items).

Conflicts of Interest
The authors declare that they have no conflicts of interest.