Abstract
Water inrush in underground mines is a major safety threat for mining personnel, and it can also cause major damage to mining equipment and result in severe production losses. Water inrush can be attributed to the coalescence of rock fractures and the formation of water channel in rock mass due to the interaction of fractures, hydraulic flow, and stress field. Hence, predicting the fracturing process is the key for investigating the water inrush mechanisms for safe mining. A new coupling method is designed in FRACOD to investigate the mechanisms of water inrush disaster (known as “Luotuoshan accident”) which occurred in China in 2010 in which 32 people died. In order to investigate the evolution processes and mechanisms of water inrush accident in Luotuoshan coal mine, this study applies the recently developed fracture-hydraulic (F-H) flow coupling function to FRACOD and focuses on the rock fracturing processes in a karst collapse column which is a geologically altered zone linking several rock strata vertically formed by the long-term dissolution of the flowing groundwater. The numerical simulation of water inrush is conducted based on the actual geological conditions of Luotuoshan mining area, and various materials with actual geological characteristics were used to simulate the rocks surrounding the coal seam. The influences of several key factors, such as in situ stresses, fractures on the formation, and development of water inrush channels, are investigated. The results indicate that the water inrush source is the Ordovician limestone aquifer, which is connected by the karst collapse column to No. 16 coal seam; the fracturing zone that led to a water inrush occurs in front of the roadway excavation face where new fractures coalesced with the main fractured zone in the karst collapse column.
1. Introduction
Water inrush has long been a hydrogeological hazard in coal resource exploitation. The hydrogeological conditions of coal mines in North China and East China are complex, and the threat of water inrush disaster is becoming more and more serious with the increase of mining depth. In the past 10 years, there have been 52 major water inrush accidents, with a direct economic loss of more than RMB 3 billion (US$46.3 million) (Guo W 2018). Water inrush has become a major hidden danger in coal mine safety production. The water inrush accident at Luotuoshan coal mine on March 1, 2010 [1], is one of the most deadly accidents in record, which resulted in a loss of 32 lives and direct economic loss of up to RMB 48.53 million (US$7.5 million). In recent years, with the increase in mining depth [2], the influence of high-pressure water on coal seam floors is becoming more severe, and the likelihood of water disasters is increasing.
The study of the water inrush mechanism and the prediction of mine water flood will contribute to coal mine safety and national production safety performance. Regardless of the type of water inrush, the main causes [3] are the large deformation and failure of the surrounding rock of the mining voids. Rock fracturing is the dominant mechanism of failure, especially in intact rocks where explicit fracturing in the zones of geological structures or breakage of the roadway floor, rather than plastically deforming intact rock, are dominant [4]. Therefore, understanding the evolution processes and mechanisms of water inrush channels, controlled by fracture propagation, could be the key to preventing inrush events.
Water inrush events have been studied by using theoretical analysis [5–7], numerical simulation using continuum approaches [8], and physical simulation [9, 10]. However, in the above studies, coal water inrush was analysed from large-scale regional predictions, whereas the evolution and mechanism of water inrush channels caused by mining and groundwater were not described directly.
FRACOD is designed to simulate fracture initiation, propagation, and coalescence in elastic rock and has been well developed. Now, it also designs complex coupling processes among rock fracturing, heat transfer, and hydraulic flow. Since its initial development in the 1990s, FRACOD has been used in many studies including the well-known ÄSPÖ Hard Rock Laboratory’s Pillar Spalling Experiments (APSE) in Sweden [11], the Decovalex International Collaboration Project [12], the Mizunami Underground Research Laboratory (MIU) Investigations in Japan [13], and the State Key Laboratory of Mining Disaster Prevention and Control in China [14]. In recent years, the focus of FRACOD development has shifted to the coupling between rock fracturing, fluid flow, and thermal processes. The coupling between rock fracturing and fluid flow was achieved using an explicit approach, which employs a time marching iteration scheme for both fluid flow and fracture propagation [15]. The new function is capable of simulating fluid flow in complex fracture networks and fracture movement and propagation driven by fluid pressure.
This paper describes a study using the numerical code FRACOD to investigate the mechanisms of the Luotuoshan water inrush disaster. The theoretical background of FRACOD and its recent developments of the fracture-hydraulic flow coupling are first introduced. Then, we present a case study to investigate the mechanism causing and evolution of the water inrush channels of the Luotuoshan accident with a focus on how the fracture initiation and propagation linked the excavation with the collapse column. The case study is conducted based on the actual geological conditions of Luotuoshan mining area, and various materials with actual geological characteristics are used to simulate the rocks surrounding the excavation of coal seam. The influences of various key factors are analysed, such as in situ stresses, fracturing on the formation, and development of water inrush channels. By introducing random small cracks in the karst collapse column [16], a new way of modeling the initiation and propagation of water inrush channels has been adopted. The coupling process of fracture-hydraulic flow during the formation of water inrush was studied, and the fracturing pattern in surrounding rock of a roadway heading was modelled considering the interaction of water pressure and rock stresses.
2. Analysis of Water Inrush in the Luotuoshan Coal Mine
2.1. The Luotuoshan Water Inrush Accident
The Luotuoshan coal mine is located in Wuhai City in the Inner Mongolia Autonomous Region of China, and the original design production capacity was 1.5 Mt annually. The mine access includes two inclined shafts and one vertical shaft (main inclined shaft, auxiliary inclined shaft, and vertical ventilation shaft) and the mine has two mining levels (+870 m, +920 m).
On March 1, 2010, the main return roadway of the No. 16 longwall panel at the +870 m level of the Luotuoshan coal mine experienced a sudden increase in water inflow. At 5:50 pm, water spraying and floor heave were first observed at the cutting face of the No. 16 main return roadway. About 90 minutes later, the mine water inrush occurred accompanied by a loud noise. The average rate of inrushing water reached 65000 m3/h, the water pressure was estimated to be 4 MPa, and the total volume reached 75000 m3. Because the quantity of maximum inrushing water was far greater than the drainage capacity of the mine, the mine flooded, and 77 people were trapped underground. The underground water level rose to +1093 m after stability, as shown in Figure 1. This water inrush accident ultimately led to a loss of 32 lives and a direct economic loss up to RMB 48.53 million (US$7.5 million).
2.2. Analysis of Mine Hydrogeology Conditions
The strata of Luotuoshan coal mine are, from bottom to top, Ordovician (), Carboniferous (86 m), Permian (103 m), Triassic (211 m), and Quaternary (6 m) systems. The main mining coal seams are No. 9 and No. 16, the average mining thickness is 4.33 m in the No. 9 seam and 5.11 m in the No. 16 seam. The groundwater regime of mining area at Luotuoshan consists of three subsystems: pore aquifer formation of Quaternary, sandstone fracture aquifer formation, and limestone karst aquifer formation of Ordovician system.
Based on pumping-out test data of limestone karst aquifer of Ordovician system in Luotuoshan [17], the aquifer thickness ranges from 20.04 to 26.62 m, with an average value of 23.33 m, and the level of groundwater ranges from +1259.21 to +1269.49 m, with an average value of +1264.35 m. The water inflow is estimated at 0.0291 to 0.0311 L/s, the unit-inflow ranges from to L/s·m, with the permeability coefficient of 1.62~2.05 mm/d, and water pressure of around 4 MPa. Thus, the aquifer has weak water yield property, ground water recharge condition is poor with incomplete development of karst fissure, and hydrogeological boundary condition is simple. The location of water inrush accident, main return roadway of No. 16 longwall panel, is 33 m from Ordovician.
2.3. Characteristics of the Collapse Column and Stratum
Formation of a karst collapse column is a dynamic geological phenomenon in which the rock mass or deposits overlying a karstified zone subside along a karst cavity, resulting in a collapse pit or sinkhole [18]. The karst collapse column has a close relationship with the extent of karst development, the character and the thickness of overburden, and the dynamic condition of underground water.
Based on the detailed geological prospecting work conducted after the water inrush accident [19], the spatial distribution of the collapse column is shown in Figure 2(a) and Figure 2(c). Z1, Z2, and Z3 are geological drilling; Z1-1, Z2-1, and Z3-1 are supplementary drilling; T1 is exploratory drilling of collapse column after the accident. The engineering geological drilling exploration, as shown by T1, is carried out in the water inrush area, and stratigraphic structure (Figure 2(b)) is revealed by borehole near the point of water inrush. The collapse column consists of two parts, the main cave and secondary cave, and the physical properties of the rocks surrounding the roadway and coal are shown in Table 1. The main cave is less than 10 m in front of the working face of the No. 16 main return roadway, and the secondary cave is buried at the bottom of the coal seam, which has well-developed fractures that allow water flow.
(a)
(b)
(c)
From the plane view and section view of the collapse column, it is believed that the water inrush channel is formed by the fractures connected with the small top section of the collapse column, based on features of water inrush and bursting mass. The source of water irruption comes from the aquifer of Ordovician.
3. Fracture-Hydraulic Flow Coupling Models in FRACOD
FRACOD is designed to simulate fracture initiation, propagation, and coalescence in elastic rock. The basic theories and various applications to engineering problems can be found in Shen et al. [20–22] and Stephansson et al. [23]. At present, most of the software used to simulate rock crack propagation adopts the methods of finite element and discrete element [24, 25], which is different from FRACOD software. FRACOD is a two-dimensional boundary element code which follows the boundary element method (BEM) principals. The code employs the displacement discontinuity method (DDM) to predict the explicit fracturing process including the possibility of fracture initiation and possibility and the path of fracture propagation.
3.1. Theoretical Background of FRACOD
Tensile and shear failure are common in rock masses. When simulating fracture process, it is needed to consider mode I (tension), mode II (shear), and mixed mode I+II propagation. G-criterion (Maximum Strain Energy Release Rate Criterion) in its original form is not directly suitable for predicting both mode I and mode II fractures. Consequently, the original G-criterion had been improved and extended by Shen and Stephansson [26] to form a new criterion, the F-criterion. The principle of the F-criterion can be stated as follows.
In an arbitrary direction () at a fracture tip, is defined by where and are the critical strain energy release rates for mode I and mode II fracture growth; and are strain energy release rates due to the potential mode I and mode II fracture growth of a unit length. These are calculated numerically with these steps: (1) compute the strain energy based on the original crack; (2) add a fictitious crack increment in direction , and compute the strain energies for the crack growth modes I and II based on this new fictitious crack; (3) divide the difference between the new and original strain energies by the incremental length to obtain the approximate values of and .
Let reach the maximum value at . The F-criterion states that if , fracture propagation will occur in the direction .
3.2. Fracture-Hydraulic Flow Coupling
It is known based on the field and laboratory observations that the interaction characteristics between pressurized water and surrounding rock have a major effect on controlling floor water inrush. In a fractured rock mass, water flow and change in fracture networks are interrelated. Rock fractures will enhance water flow by creating new flow channels, whereas water pressure may stimulate fracture growth. In order to study the water inrush issue in Luotuoshan mine, understanding and simulating the coupled processes between fracture mechanical response and water flow are necessary.
There are two fundamental approaches in modeling the coupled fracture-hydraulic flow (F-H) processes in fractured rock. One is the implicit approach and the other is the explicit. In the first approach, related fluid flow equations [27], such as Darcy’s law, are solved together with mechanical equations for rock medium. In the second approach, however, both fluid flow and mechanical response are simulated by using a time marching iteration process. There are several well-known commercial codes, such as UDEC [28], that are based on the explicit process.
In order to simulate the coupled processes, a new coupling method, which is an explicit approach, is designed in FRACOD, and it uses the displacement discontinuity (DD) method with an iteration process to simulate rock deformation and fracture propagation. In parallel, it also uses the cubic law combined with an iterative process to model the fluid flow processes in fractures. Comparing with the implicit approach, both the mechanical and fluid flow calculations with this explicit approach are mathematically simple to adopt complicated boundary conditions and changing model conditions. But it often needs longer computational time for flow solution.
3.2.1. Numerical Considerations
In general, there are two ways of water inrush in mining face: water not only flows in the fracture channels but also leakages from the fracture channels to the surrounding intact rock. The paper considers both the flow in fractures and the leakage in intact rock, as shown in Figure 3. The path of water inrush is dominantly the fracture channels, the so-called “channels” water inrush or fracture water inrush, and in the rock it is seepage water inrush.
During the mechanical simulation using DD method, a fracture is discredited into a number of DD elements, as shown in Figure 4. In the flow calculation, each DD element is considered as a hydraulic domain, and the adjacent domains are connected hydraulically. Fluid may flow from one domain to another depending on the pressure difference between the two domains.
3.2.2. Iteration Scheme
A coupled F-H process is illustrated in Figure 5, and the iteration steps are described below.
Step 1. Water flow occurs between fracture domains, and water leaks into the rock matrix. The water flow between the fracture domains is calculated using the cubic law, that is, the flow rate between two adjacent domains is calculated using
where is the fracture aperture; is the domain (element) length; is the fluid pressure difference between the two domains; and is the dynamic fluid viscosity.
The leakage between a fracture domain and the rock matrix is described by porous flow, but the porous flow process is often complex and it depends on the varying boundary conditions, initial fluid pressure, and flow history [29]. To simplify the simulation, it is assumed that the fracture is effectively enclosed by a fictitious layer of intact rock of limited thickness at each side and a constant fluid pressure exists on the outer boundary of the rock layer. The flow leakage between a fracture domain and the surrounding intact rock is estimated using
where is the permeability of rock; is the effective leakage distance; is the domain fluid pressure; and is the initial pore pressure minus the pressure at the outer boundary of the simulated rock layer.
Step 2. Water flow causes pressure changes in the water domains. After a short time duration , the domain pressure will change as the water flowing in or out of a domain changes the fluid volume. The new domain pressure is calculated using where is the fluid bulk modulus; is the domain volume; and is the time step.
Step 3. Changes in water pressure cause fracture deformation. The fracture domains are calculated using the DD method, where the water pressures in the fracture domains are the input boundary stresses. After considering the water pressure in the fracture domains, the system of equations for calculating the element displacement discontinuities is given as where and are the in situ stresses in the tangential and normal directions of the fracture surface; and are the displacement discontinuities in the shear and normal directions; , , , and are the influence coefficients of the displacement discontinuity on the tractions on the fracture elements; and and are the shear and normal stiffnesses of the fracture.
Step 4. The fracture deformation changes the domain volume, thus changing the water pressure in the domain. The new domain pressure is calculated using where is the initial domain pressure before fracture deformation; is the new domain pressure after fracture deformation; and is the change in fracture aperture.
The new domain fluid pressures are then used to calculate the flow rate between domains in Step 1. Steps 1 to 4 are repeated until the desired fluid time is reached and a stable solution is achieved.
3.3. Validation and Application of FRACOD
Over many years, FRACOD has been used to model the rock engineering problems and has been proven to be useful in predicting borehole breakouts [30], stability of large shaft and galleries [31], pillar spalling [32], rock mass permeability change due to fracturing, and fundamental creep behavior of rock samples [33].
To demonstrate the fracture propagation mechanisms discussed in the coupled F-H code of FRACOD, a single fracture fluid flow example case has been modelled. This case only tests the fluid flow function with mechanical coupling; hence, the single fracture linking the two holes has a constant initial aperture of 10 μm. Two boreholes are of a diameter of 1 m and are 10 m apart. The left hole is the inlet hole and has a constant fluid pressure of 6 MPa; the right hole is the outlet hole and the fluid pressure is 0, as shown in Figure 6. It was excavated at a depth of about 100 m below the ground surface, with horizontal (Syy) and vertical (Sxx) in situ stress and of about 2.5 MPa. The single fracture was set as the monitoring line to monitor the change characteristics of fluid pressure and flow rate at different times. The key mechanical and fracture properties used in the simulation test are shown in Table 2.
The modelled fluid pressure distribution with time in the single fracture is shown in Figure 7(a). The results demonstrate a dynamic process of fluid flow from the inlet hole to the outlet hole in 0.2 s, 0.5 s, 1.5 s, and 10 s. Affected by the water pressure between two holes, the fluid flows from the inlet to the outlet. After 0.2 s, 0.5 s, and 1.5 s, the pressure distribution is still nonlinear before 1.5 s. After about 10 s, the fluid flow in the single fracture has reached a steady state solution and the pressure is linearly distributed. The fluid flow distribution is shown in Figure 7(b). After about 10 s, the fluid flow has reached a steady state and the flow is linearly distributed. Before reaching steady stage, the flow near the inlet is larger than that in other areas, and the flow gradually decreases with the distance increase from the inlet. The results of the tests show that the F-H code of FRACOD is working properly according to the expectation.
(a)
(b)
We use a simple model for the main return roadway which is located in the No. 16 coal seam (Figure 8). The simulated heading face of the main return roadway is arranged at the top of the No. 16 coal seam; the width and height are 5.22 m and 3.66 m, respectively, and it was excavated at a depth of about 414 m below the ground surface, with vertical and horizontal in situ stress of about 10.3 MPa. In this study, the mechanical properties of the simulated rocks are the same as those in the field. Based on the physical properties of coal rock, Young’s modulus is 10 GPa, tensile strength is 1.0 MPa, Poisson ratio is 0.25, the angle of internal friction is 35°, and cohesion is 2 MPa. In this case, fracture initial and residual aperture is assumed to be 10 μm, and fracture shear and normal stiffness is assumed to be100 GPa/m.
(a)
(b)
(c)
(d)
The modeling was conducted using symmetry fractures to reduce the calculation time. In this case, we can clearly see that this approach has proved very effective in simulating the behavior of fracture initiation, propagation, and coalescence in surrounding rock of the main return roadway. Because the natural stress field is redistributed due to the excavation, stress concentration areas appear in rock surrounding the roadway that causes roadway deforming or destructing. Small fractures were first observed in the roadway’s sidewalls near the corners due to tensile stress and low confinement (Figure 8(a)). High stress concentration occurred at the roof and bottom corners; the maximum value of the major principal stress is about 45 MPa (Figure 8(c)); and the initial fractures started to propagate in these areas and extended in sidewalls, roof, and floor (Figure 8(b)). Although the propagation of fractures was dominated by tensile failure, shear failure played an important role in the path of propagation. The maximum value of roadway displacement occurs at the surface of the roof (Figure 8(d)). This modeling case demonstrates that the code is working properly according to the expectation.
3.4. FRACOD Models to Simulate the Water Inrush
There are three factors affecting mine water inrush: water source, water inrush channels, and water inrush condition. Among them, the presence of a water inrush channel is a key factor affecting mine water inrush. Previous studies mostly concentrated on the floor damage characteristics and did not provide a detailed description of the formation of water inrush channels under the joint action of the mining process, geological structure, and water pressure. The aim of this article is to elucidate the formation of a water inrush channel during the process of tunnel excavation.
The latest version of FRACOD with M-H coupling functions is used to simulate the water inrush in this paper. A roadway excavation model was established along the direction of face advance (see Figure 9), and the spatial parameters of each structure are the same as those of the actual geological data in this model. The simulated main return roadway in the No. 16 seam is 36 m long and 3.66 m high, excavated at a depth of about 414 m below the ground surface. The roadway was located in the upper part of the No. 16 coal seam, and the roadway floor is 38 m above the limestone karst aquifer. The main cave of collapse column is in front of roadway working face, with the distance less than 10 m, and the secondary cave is buried in the bottom of the coal seam.
Note that the modeling was done in two-dimensional (2D). The 2D representation of a roadway along its axis direction is problematic because it assumes the infinite width of the roadway, which is clearly an unrealistic assumption. Nevertheless, the 2D model is considered to be more conservative than reality because it can lead to a higher stress concentration in front of the roadway face. The mechanisms of rock fracturing in the 2D model should still be reasonably representative of reality.
During this numerical study, , the model includes two parts: surrounding rock and collapse column, as shown in Figure 9. The surrounding rock mass is treated as strata of different homogeneous, isotropic, and elastic media, and geostatic stress is considered. The rock mass is assumed to be impermeable and nonporous due to the code limitation in FRACOD. It needs to emphasize that the collapse column should not be directly simulated as one continuous geological body, because it requires a much more complex model and significantly longer computational time to achieve convergence for curvilinear boundary of collapse column. In part of the collapse column, these cracks with random orientations are assumed to have a hydraulic aperture that enables them to hold water pressure. However, due to the limitation in modeling capacity, it is impossible to define the boundary of the collapse column and use the actual sizes of the microcracks. Hence, as a simplification, we used cracks with a length of 1.4 m and a crack density of 0.5 cracks per square meter. The fractures near the boundary of the collapse column are defined at the end point of the fracture line by a tip element, and the fracture initiation (element size of 0.5 m) is allowed in the internal rock mass and at the strata boundary.
The mechanical properties of the simulated rocks and coal are given in Table 3. In this study, Young’s moduli represent the lower end values of a range of measurement results using different methods. This treatment is designed to consider the scale effect from laboratory scale to field scale, because the rock mass’ Young’s modulus tends to be significantly lower than that of a small rock specimen in the laboratory scale.
4. Modeling Results and Discussions
4.1. A Description of the Water Inrush Event That Occurred
Because the natural stress field changes due to roadway excavation, an area of stress concentration forms in the rock surrounding the roadway. Meanwhile, the fracture and damage of the rock mass significantly change the permeability of the surrounding rock, thus causing water inrush accidents in coal mines. In this paper, the factor of safety (FoS, see [33]) is used to analyse failure pattern of surrounding rock at the main return roadway of the No. 16 seam. FoS is a term describing the load carrying capacity of a system beyond the expected or actual loads. Essentially, the factor of safety measures how much stronger the system is than an intended load.
In FRACOD, the total FoS has two main components: shear FoS and tensile FoS. The total FoS is a ratio of maximum shear/tensile strength of materials to intended load with the actual shear/tensile values of materials. If the value is greater than 1, surrounding rock could take additional load before it fails by this definition. If the value is less than 1, the surrounding rock will fail, and deformation in the rock mass will produce fracture propagation, allowing water inrush.
Figure 10 shows the total FoS distribution in the surrounding rock mass at different times. The failure zone of surrounding rock is gradually increasing with time near the roadway face, especially between the main cave of collapse column and the roadway face. After cycle 9 (Figure 10(a)), the failure zone between the collapse column and the roadway face reaches the main cave, and the failure zone in the roadway floor started to expand towards the secondary cave. After cycle 12 (Figure 10(b)), fluid flow in the rock mass in the failure zone of the main cave may cause localized tensile stress, which further drives the preexisting fractures to propagate. After cycle 14 (Figure 10(c)), the water pressure drove the short fractures in the main cave to propagate towards the roadway face due to expansion of the failure zone and coalescence with the main fracture. At cycle 24 (Figure 10(d)), the major fractures coalesce, water inrush channels form, and water rushes into the roadway. It can be clearly seen that the water inrush location is in the area at the front of the roadway heading, as the fracture channels coalesced with the main cave in the collapse column.
(a)
(b)
(c)
(d)
4.2. Stress of the Surrounding Rock at the Roadway
Stress distribution characteristics in rock surrounding the main return roadway of the No. 16 depend on many factors, such as rock mechanical properties (e.g., Young’s modulus) and the effect of geological structures. The normal stresses of the surrounding rock are shown in Figure 11.
(a)
(b)
As shown in Figure 11(a), the simulated normal stresses (Sxx) along monitoring lines Lv1-4 (Figure 4) for cycle 14 follow a similar trend and have very close values at the lower part of the lines but large differences at the upper part of the lines in the region in front of the roadway face. This is because the rock in the region at the lower part of the monitoring lines is relatively intact, while fractures are well developed in the region of the coal seam at the front of the roadway face (see Figure 10(c)). The maximum normal stress is observed at the top end of monitoring line Lv1, near the corner of the roof. At the upper part of the monitoring lines, the normal stress decreases with increasing distance from the roadway face. The variation in the normal stress along the upper part of line Lv1 is different from the variations in the other three lines. The variation in stress far from the roadway area influenced by mining always exhibited a hysteresis, unlike that of the rock mass near the same area. Hence, the distribution within the high-stress zone was not uniform in the pillar because of fracture initiation and propagation. The zone was first produced in front of the roadway face, and then, fractures constantly initiated over time. As a result, more water could flow into this area due to the greater number of fractures.
Geometrically, the abovementioned part is similar to the top part of a vertical monitoring line: both are near the corners. The maximum normal stress (Syy) is at the end near the floor corner of the monitoring line Lh1 (Figure 10(b)). The variations in the normal stress in the region at the right end of the monitoring lines are more complex. This may be due to the effect of the secondary cave of the collapse column (Figure 9). Some complex variations can also be seen in this region in Figure 10(c) and Figure 12(c). Nevertheless, the predicted results in Figure 9 indicate the potential impact of fracture initiation, propagation, and coalescence.
(a)
(b)
(c)
(d)
4.3. Fracture Initiation and Propagation
4.3.1. Fracturing in the Coal Pillar between the Collapse Column and Roadway Face
The FRACOD model predicted that a major fracture initiated at the top corner of the roadway face and propagated towards the collapse column. Meanwhile, the initial short fractures at the collapse column started to propagate towards the roadway face. Although the number of those fractures was very small, they play an important role in the formation of water inrush channel due to the fracture coalescence. The fracture propagation was caused because the high stresses in this area exceeded the strength of the coal and rock. The high stresses are due to the fluid pressure and stress redistribution from the excavation of the roadway. The propagation of these fractures was dominated by shear failure in this case, however, tensile failure played an important role in the path of propagation.
In the early stages of modeling, it was observed that the new fractures in the main cave initiated in a direction approximately parallel to the preexisting fractures (Figure 12(a)). As time progressed, the water pressure drove the short preexisting factures in the main cave to propagate towards the roadway face (Figure 12(b)), and eventually, these fractures coalesced with the factures initiated in the coal pillar and those from the roadway face (Figure 12(c)). The newly linked fractures then propagate under high fluid pressure, eventually forming several long-inclined fractures and then water inrush channels (Figure 12(d)).
4.3.2. Fracture in Coal Floor
The fracture development in the coal floor has two effects: one is the broken floor strata, and the other is the propagation of the preexisting fractures of the secondary cave in the subvertical direction. However, these two fracture effects did not interact during the whole course of modeling.
Based on simulation results, it is judged that the source of the water inrush accident in Luotuoshan coal mine is the Ordovician limestone aquifer, which is connected by the collapse column to the No. 16 coal seam; water inrush occurred in front of roadway heading because fractures in the main cave in collapse column driven by excavation stress and high-water pressure propagated and coalesced with those initiated at the roadway corner.
4.4. Formation and Development of Water Inrush Channels
The above simulation results indicate that water inrush occurs in front of the roadway excavation face where new fractures coalesce with the main fracture zone in the karst collapse column. A new numerical simulation of a roadway is conducted to study the formation and development of water inrush channels before and after fracture coalescence. The mechanical properties and burial conditions of the simulated rocks are the same as the verification model. In this case, the model consists of the cross section of the roadway and some preexisting fractures at one side, representing the collapse column (see Figure 13).
(a)
(b)
Fractures filled with high-pressure water were observed to propagate in the collapse column towards the roadway face (Figure 13(a)), and a major fracture without water at the top corner of the roadway face propagated towards the direction of the collapse column. The blue arrows on the preexisting fractures represent the relative magnitudes of the fluid pressure. After the fractures coalesced (Figure 13(b)), water inrush channels formed, and high-pressure water gushed into the roadway face. The high pressure from the water released, and numerically, the water pressure decreased to zero. The newly initiated fractures are not predicted to propagate under the stress conditions, and no water inrush channel has been predicted to form. At this stage and after, no stress-induced fracture initiations are predicted because the stresses are released rather than elevated due to excavation, and the surrounding rock mass of the main return roadway of No. 16 generally stabilizes.
5. Conclusions
Understanding the coupled fracturing processes from the underground water is important for preventing coal water inrush. To investigate the evolution process and mechanisms of the major water inrush accident in Luotuoshan coal mine, this study applies the recently developed fracture-hydraulic flow coupling function of FRACOD. The F-H coupling was achieved using an explicit approach, and it employs an iteration scheme for both fluid flows and fracture propagation. The new function is capable of simulating fluid flows in complex fracture networks and fracture movement and propagation driven by fluid pressure.
The case study conducted for the Luotuoshan accident using a FRACOD model indicates that the water pressure drove the short factures in the main cave to propagate towards the roadway face and eventually coalesce with the factures in the coal pillar. These fractures propagated and coalesced with other newly initiated short fractures and those initiated from the roadway corner, eventually forming several long fractures and hence a water inrush channel. Based on the simulation results, it is determined that the water source of the accident was the Ordovician limestone aquifer beneath the mining seam, which is connected by the collapse column to the No. 16 coal seam. Water inrush occurred when the roadway face was at some distance from the main cave in the collapse column, where fracture initiation and propagation finally linked the roadway to the main cave and formed water inrush channels.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declared that they have no conflicts of interest to this work.
Acknowledgments
This research was financially supported by the National Natural Science Foundation of China (51974173), Key Research and Development Plan of Shandong Province (2019GSF111024), Open Fund of State Key Laboratory of Water Resource Protection and Utilization in Coal Mining (SHJT-17-42.14), and Shandong Province’s Taishan Scholar Talent Team Support Plan for Advantaged & Unique Discipline Areas. The authors also thank Dr. Jingyu Shi from CSIRO Energy for his help in polishing the wording and writing.