Abstract
Although the mechanism and influence of fault water inrush have been widely studied, there are still few studies on the migration of filling particles and the evolution process of seepage characteristics within faults. In this work, the coupling effects of water flow, particle migration, and permeability evolution are considered synthetically, and the evolution model of seepage characteristics with multifield coupling is established. This model was used to investigate the evolution process of water inrush within faults and the effects of water pressure, initial effective porosity, and initial permeability on water flow rate. The results show that the evolution of seepage characteristics can be divided into three phases: (i) low velocity seepage, (ii) drastic changes with substantial particle migration, and (iii) steady-state water flow. The multifield coupling causes the effective porosity, permeability, flow velocity, and particle concentration to accelerate each other during the dramatic phase. Moreover, the increases in initial water pressure, initial porosity, and initial permeability have different degrees of promotion on the water flow rate. Finally, the simulation results are approximately the same as the data of water inrush in the mining area, which verifies the correctness of the evolution model established in this work. This work provides new approaches to the evolution process and prevention of water inrush in faults.
1. Introduction
Affected by complex geological and hydrological conditions, water inrush disaster has been a serious challenge in underground engineering such as coal mines in China [1–3]. The existence of the geological structure usually functions as an underground water flow path, thus posing a great threat to the mining production safety [4, 5]. Fault structure is the main factor leading to water inrush in the mine pits. According to statistics, 80% of water inrush accidents are related to faults [6]. For instance, the 9101 working face in the Yangzhuang coal mine of the Feicheng Mining Bureau was between two walking faults, and the fault structures communicated the mining floor with the Xujiazhuang limestone aquifer and Ordovician limestone aquifer, which led to the occurrence of water inrush in the working face two times. Among them, the water pressure of the second water inrush reached 8.5 MPa, with a water flow rate of nearly 600 m3/h. And the water flow rate sharply increased to 4409 m3/h after 16 hours, which has greatly exceeded the drainage capacity. After 19 hours, the depth of water in the shaft station reached 0.7 m. Finally, the mine was all submerged due to the lack of disaster resistance [7]. Thus, fault structure, as an important water inrush channel, is of great significance to the prevention and control of fault water inrush hazards in coal mines.
Aiming at the mechanism of water inrush caused by geological structures such as faults, the scholars have conducted a large number of experimental investigations. Through laboratory tests, Wu et al. [8] systematically analyzed the deformation and destruction patterns of fault materials under different moisture content, confining pressure, and loading modes. It was pointed out that the intensity of fault materials will be enhanced with the increase in confining pressure, but reduced with the increase in moisture content, and the deformation and destruction characteristics of the fault materials have a close relation with its materials and loading time. Li et al. [9, 10] developed a large-scale testing system for the measurement of coupled seepage and triaxial stress. The seepage failure model of the filling medium was further revealed. The catastrophe evolution mechanism of water inrush for the filling fault had been researched, and the variation characteristics of fault occurrence, filling ability, and water conductivity have been analyzed. Zhang et al. [11] simulated physically the entire process of crack formation, concealed fault propagation, and evolution of a water inrush channel with high-pressure water directly beneath the mine floor, and the results indicate that water channels are mainly caused by the connection between tectonic rock zones and coal floor cracks, which are the direct cause for water inrush. Miao et al. [12, 13] carried out a series of tests on the permeability characteristics of the fractured rock mass and obtained the variation rules of the permeability of the fractured rock mass under different lithology and stress states. Ma et al. [14–17] carried out a series of seepage experiments on granular rock. The effects of the original particle size distribution (PSD) on the hydraulic properties were investigated, and the evolution of porosity and permeability was calculated and verified in detail.
Many scholars also investigated the theoretical analysis and prediction of fault water inrush. Using digital image processing technology, Zhao et al. [18] developed MATLAB function for threshold segmentation of the section image of broken rock mass and established a three-dimensional digital analysis model. Based on these models, interconnected networks of spatial structures and spatial distribution features of stress, seepage, and speed are identified for cataclastic rock masses. Liu et al. [19] established a coupled thermal-nonlinear-hydraulic-mechanical (THM) model for fault water inrush to study the water-rock-temperature interactions and predict the fault water inrush. Yao et al. [20] presented an inverse velocity theory to predict the water inrush time under different geological and flow conditions. The results showed that the inverse velocity theory was capable of predicting the occurrences of water inrush under different conditions, and the time of water inrush had a power relationship with the rock heterogeneity, water pressure, and initial particle concentration. Song et al. [21] analyzed strata’s movement and ground pressure caused by mining’s influence on karst water pressure and fault pillar. It was pointed that karst water pressure’s increase has an important relationship with overburden strata’s movement caused by mining. Especially for the closed karst cave, water pressure will leap under the effect of abutment pressure. In order to investigate the activation failure characteristics of the mining floor fault above a confined aquifer and predict the water inrush from the activated floor fault, Sun and Wang [22] and Jian et al. [23] used microseismic monitoring technique to monitor the mining floor fault.
The above results are of great significance to reveal the mechanism and prediction of fault water inrush. However, the seepage characteristics of faults are the determinants to control the water conductivity of faults and the degree of damage. Permeability is a direct reflection of the water conductivity of faults. In the process of seepage, the filling particles are constantly migrated and lost due to the erosion of water flow, resulting in changes in porosity and permeability as well as particle concentration. However, the changes in overall permeability will react on the seepage field. It is a process of mutual influence and promotion. Therefore, the coupling effect will be an important reason to reveal the process of fault water inrush. The above correlation analysis rarely involves the influence of the coupling effect between particle migration and seepage characteristics, and the evolutionary process of the interaction has not been systematically studied.
In this work, the author comprehensively considered the coupling effect of water seepage and particle migration in the process of fault seepage and established a multifield coupling model for the evolution of the seepage characteristics. Combined with the geological conditions in the mining area, the governing equations will then be implemented into COMSOL Multiphysics software. The heterogeneity of rocks in faults will be determined using the Weibull distribution, and the evolutionary feature of seepage characteristics, including effective porosity, permeability, particle concentration, and seepage velocity, will be obtained accordingly. The process of seepage characteristics interaction will be analyzed, and the variation of flow rate under different natural conditions will be discussed in detail. The evolution process of the seepage characteristics will be simulated and investigated, and the mechanism of water inrush from the perspective of multiple coupled system movement will be revealed.
2. An Evolution Model of Seepage Characteristics
2.1. Hydrogeological Structural Features of Faults
The hydrogeological structure of the fault is shown in Figure 1, which can be regarded as composed of broken solid medium (skeleton and filling material), liquid medium (fluid) in the pores and fractures, and fine filling particles in the liquid medium. The surrounding rock is relatively hard, and the bottom is confined aquifer. When mining reveals broken geological structures such as water-carrying faults, seepage erodes the filling particles in the framework of broken rock mass and causes migration as well as changes in the permeability characteristics of faults, which may lead to water inrush hazards.
2.2. Governing Equations
2.2.1. Assumptions
In order to establish the evolution model of seepage characteristics and investigate the mechanism of fault water inrush in coal mines, the following assumptions referring Yao et al. [20] have been made: (1)The fluid and particles in faults are incompressible(2)The suspended particles share the velocity field with the fluid(3)The effect of erosion on the permeability change is in proportion to the change of particle concentration in the fluid(4)Fault fillings can be regarded as porous materials
2.2.2. Mass Conservation Equations
A representative element for a typical fault-filling microstructure is shown in Figure 2. Suppose that the coordinate of the element is , then its volume (m3) can be expressed as follows:
According to the above analysis, the element is composed of three kinds of media: (1) solid medium with a volume of (m3), (2) fluid medium with a volume of (m3), and (3) particles suspended in the fluid with a density of (kg/m3), a volume of (m3), and a mass of (kg), where the subscripts “s, f, and fs” represent solid medium, fluid medium, and particles in the fluid, respectively. Then, the pore volume (m3) and porosity (%) of the element can be expressed as follows [24]:
The volume concentration (%) and mass concentration (kg/m3) of the suspended particles can be expressed as follows:
The migration of particles in volumetric elements can be regarded as the result of the combined action of advection-diffusion. As particles with extremely small size migrate in high-speed flow fields, the tiny diffusion effect can be neglected. Along the -axis, the average velocity of particles can be defined as (m/s), and the particle mass in the inflow of the element caused by advection in the -direction per unit time can be expressed as follows:
The particle mass in the direction per unit time can be expressed as follows:
For the unit element, the total particle mass per unit time can be summarized as follows:
The decrease in element mass per unit time can be expressed as follows: where (kg/m3/s) is the particle mass that migrates into the fluid from a unit element in a unit time, and is written as follows:
According to the mass conservation law, by combining equations (7) and (8), it can be summarized as follows:
Assuming that the solid particles are incompressible, and the volume of the element body does not change with time, combining with equations (3) and (4), equation (10) can be summarized as follows:
Equation (11) can be expressed as a vector expression as follows: where .
Substituting equation (9) into equation (12) yields
Similarly, the mass conservation equation of the fluid can be obtained as follows:
2.2.3. Water Seepage
Darcy’s law describes linear seepage in porous media, and it is derived from balance of momentum for the fluid phase and from a constitutive equation for the fluid-solid interaction force (i.e., the seepage force) [25]. In general, the transient form of Darcy’s law can be expressed as follows: where (kg/m3) is the density of the fluid, and (kg/(m3·s)) illustrates the mass sources of the element. The equation for Darcy velocity can be expressed as follows [24]: where (m/s) is the Darcy velocity, (m2) denotes the permeability of the element, (Pa·s) is the dynamic viscosity of the fluid, (Pa) denotes the water pressure, and is the unit vector in the direction of gravity.
2.2.4. Porosity Evolution Equation
Porosity, as an important component of seepage characteristics, affects the evolution of seepage channels and permeability. Sakthivadivel and Irmay [26] investigated the erosion effect for porous media by using both experimental and theoretical methods. Vardoulakis et al. [27] summarized the following: the evolution of permeability was affected by the porosity and particle concentration as well as the seepage velocity. The porosity evolution equation of porous media is given by considering particle migration as follows [27]: where (m-1) denotes the erosion coefficient, (%) illustrates the maximum value that porosity can reach under fluid erosion, and is the absolute value of the seepage velocity (m/s). represents the velocity component in the direction of .
2.2.5. Auxiliary Equation
The velocity of the particles can be expressed as the difference between the velocity of the fluid and the velocity of the free subsidence in the static fluid as follows [24]: where (m) denotes the diameter of particles and is the dimensionless drag coefficient, which can be expressed as follows [28]: where is the Reynolds number of particles, which can be expressed as follows:
The change of permeability in porous media can be deduced according to the change of porosity, which can be expressed as follows [28]: where (m2) denotes the initial permeability, and (%) illustrates the initial porosity.
Equations (13), (14), (16), (17), and (18) together with (22) compose the coupled processes of water transport under erosion effects, as illustrated in Figure 3. To investigate the mechanism of water inrush in faults, the above governing equations will be implemented in COMSOL Multiphysics to further assess the risk of water inrush and predict the time of water inrush under different flow conditions. There are six unknowns with equations, so the set of governing equations can provide a unique solution.
3. Numerical Model Implementation
3.1. Background of the Coal Mine
In this section, a case study was conducted based on the hydrogeological conditions of a coal mine, which is affiliated with the Jining Mining Industry Group in China. The mining operation is current extracting number 3# coal seam at approximately 900 m depth. To date, a group of faults, most of which are normal faults, were found to cross the working face obliquely. The average thickness of the fault group is 80 m, and the average dip angle is about 60°. The working face and fault group position of this mine are shown in Figure 4. The working face passes through the fault group smoothly in the mining process, and the water flow rate was about 35 m3/h. However, when the working face was advanced to 560 m, the flow rate increased rapidly from 120 m3/h. By the time the working face advanced 640 m, the flow rate had increased to 380 m3/h, and the peak flow rate had reached 432 m3/h. This process shows obvious lag characteristics, which was considered as a typical coal mine fault lag water inrush.
(a)
(b)
3.2. Material Properties of Fault Fillings
Groundwater seepage has a great influence on the mechanical properties of fault zone filling, and it is the key factor that induces rock infiltration and water inrush in the fault zone. In order to obtain the lithological parameters of fault fillings, the conventional mechanical tests were carried out. The test samples were taken from the fault of the working face, collected immediately after underground exposure, and sealed on site and transported to the laboratory. The lithological parameters obtained are shown in Table 1.
In addition, the permeability parameters of the fault filling samples were determined. In this work, the hydraulic conductivity of the specimen is determined by means of a variable head test. The schematic and object of the measuring device are shown in Figure 5. The calculation principle of this measuring device can be expressed as follows [29, 30]: where (cm·s-1) denotes the hydraulic conductivity, (cm2) is the sectional area of the head tube, (cm2) illustrates the sectional area of the specimen, (s) is the measuring period, and and represent the initial and final head differences, respectively.
(a)
(b)
The seepage properties obtained are shown in Table 2. The relation between hydraulic conductivity () and permeability () can be expressed as follows:
Digital image technology provides a method for obtaining rock parameters with heterogeneous distribution [18, 31]. The effective porosity of the specimens is determined by digital image processing technology. In particular, the relationship between the pore pixel and the total image pixel is obtained through the threshold value on the slice images of the specimens, and the image processing for specimens is shown in Figure 6. The effective porosity obtained is shown in Table 2.
3.3. Numerical Model Setup
In this paper, the evolution of the seepage characteristics and the water inrush from the bottom to the top are mainly investigated, so the model is simplified into a 2D fault section (i.e., a parallelogram). According to the lithology of the mine, the bottom length and the height are 18 m and 50 m, respectively (as shown in Figure 4). The inlet water pressure at the bottom of the model was 2.0 MPa, and the outlet pressure at the top of the model was 0.1 MPa; the initial water pressure in the model was 0.1 MPa, the initial particle volume concentration was 0.01, and the initial effective porosity was 10%. Other used parameters are listed in Table 3.
Heterogeneity is a natural feature of rock. There are currently various approaches to obtaining the characteristics of the heterogeneous distribution of rock materials. For example, digital core techniques and mathematical statistics methods have both been utilized. Previous studies have shown that the heterogeneity of rock can be described by the Weibull distribution [31–34]. Thus, the Weibull distribution was selected in this study due to its effectiveness and great simplicity to obtain the heterogeneous distribution of effective porosity in the fault. The distribution probability density equation was as follows: where indicates the effective porosity, denotes the initial effective porosity, and is the uniformity index. The larger represents, the higher the level of uniformity. The initial effective porosity distribution obtained by the numerical generation method is shown in Figure 7.
4. Results and Discussion
In order to analyze the characteristics of seepage at different times, in this study, four different times were selected (i.e., 2, 4, 6, and 8 h). The specific results will be discussed below.
4.1. The Evolution of Seepage Characteristics
Figure 8 illustrates the evolution of effective porosity at different times. Figure 9 illustrates the evolution of permeability at different times. Figure 10 shows the flow velocity at different times at the outlet boundary along the -axis as marked in Figure 4(b), and the average porosity-time curve and average particle concentration-time curve are plotted in Figure 11. The specific results will be discussed below.
(a) 2 h
(b) 4 h
(c) 6 h
(d) 8 h
(a) 2 h
(b) 4 h
(c) 6 h
(d) 8 h
4.1.1. Distribution of the Effective Porosity and Permeability
From equation (22), effective porosity and permeability show a similar evolution process in Figures 8 and 9. At the initial phase, effective porosity and permeability are subject to random distribution with minor difference. As the infiltration continues, the average effective porosity increases by approximately 200%, from 10% to nearly 30% (Figure 11). However, this increase is not uniform. The variations of those with larger initial effective porosity and permeability are significant, while those with smaller initial effective porosity and permeability are slightly increased. As the distributions with large initial effective porosity and permeability preferentially increased and interconnect, some major seepage channels are formed under particle migration effect, as illustrated in Figures 8(d) and 9(d).
4.1.2. Evolution of Flow Velocity at the Outlet Boundary
As shown in Figure 10, the distribution of the flow velocity is significantly related to time. The overall increase in flow velocity is small during the initial phase, with a sharp increase in subsequent phases and a tiny increase in the final phase. The maximum hourly increase in flow velocity is approximately 1 m/s from 1 h to 2 h, 2 m/s from 2 h to 4 h, 4 m/s from 4 h to 6 h, 2.5 m/s from 6 h to 8 h, and 0.5 m/s from 8 h to 10 h. The passage of time represents the evolution of the water inrush. Therefore, the evolution process can be roughly divided into three phases according to the degree of increase in the flow rate, including (i) low velocity seepage, (ii) drastic increase in flow velocity, and (iii) steady-state water flow. In addition, the figure shows that the increases with high initial flow velocity at the outlet boundary are larger, while those with low initial flow velocity are smaller.
4.1.3. Changes in Average Effective Porosity and Average Particle Concentration
It can be seen from Figure 11 that the process of changing the average particle concentration over time has experienced a slow first phase, a drastic second phase, and a smooth third phase. Combined with the previous analysis, it can be considered that due to the small initial percolation velocity, the erosion of the water flow is weak, resulting in a slow increase in the average particle concentration; then, as the water flow velocity increases, the solid particles migrate into the fluid due to the erosion, resulting in a significant increase in the average particle concentration, and at the same time the effective porosity will also increase significantly. Though the final phase of the seepage channel and water flow rate is stable, the average particle concentration will gradually stabilize.
To sum up, the seepage characteristics have undergone a complete evolution process under the coupling of multiple fields. In the initial phase, this coupling is gradually enhanced although it is weak. As the flow rate gradually increases, the water flow continuously erodes the particles of the granular rock and migrates them. The particle concentration increases rapidly while the effective porosity and permeability also evolve. The change in permeability will in turn increase the seepage velocity and enhance the erosion, and the coupling effect causes drastic changes with substantial particle migration in the second phase. Eventually, as the porosity evolves to its maximum, the percolation channels and flow rates tend to stabilize, and the coupling is no longer significant.
4.2. Comparison between Numerical Simulation Results and Actual Data
Figure 12 plots the curves of simulated fault flow rate and actual data, and it can be seen that the variation trend of the two is basically consistent, which proves the correctness of the seepage characteristic evolution model of fault water inrush established in this paper. In connection with the above discussion, it can be concluded that the process of fault water inrush can be divided into three phases, including (i) low velocity seepage, (ii) drastic changes with substantial particle migration, and (iii) steady-state water flow. According to the actual situation, after 5 hours, the water flow will increase rapidly and water inrush hazards will occur. This time between the first phase and the second phase is called “the critical time of water inrush.” If measures such as grouting reinforcement are taken on the fault before this time, water inrush hazards can be effectively prevented. Otherwise, once the time is missed, the water flow will increase rapidly. At that point, water flow will be difficult to control and increase the risk of water inrush hazards.
4.3. Analysis of Effect Factors
According to the above analysis, the critical time of water inrush is the key time for taking measures to ensure safety. In order to study the factors affecting the evolution process, the flow conditions of water inrush were discussed. Through numerical simulation of different flow conditions, it is expected to obtain the change rules of flow rate and critical time, so as to provide a theoretical basis for water inrush prediction.
4.3.1. Effect of Water Pressure
Water pressure is the driving force of water inrush, and different water pressures must have an important influence on seepage evolution. Numerical simulation was carried out with 1 MPa, 2 MPa, 3 MPa, 4 MPa, and 5 MPa as inlet pressures. The curve of flow rate changing with time under different water pressure conditions is shown in Figure 13.
As is shown in Figure 13, water flow rate is significantly affected under different water pressure conditions. As the pressure increases, the changes in the flow rates accelerate. For example, for , the flow rate was only approximately 40 m3/h when , and it is increased to 280 m3/h when . Moreover, the critical time had a noticeable change. It can be seen that the larger the water pressure, the earlier the critical time occurs. As the pressure decreases, the interval between critical times increases. For example, the interval of critical times between and was 1.4 hours, and the interval of critical times between and increased to 2 hours. In other words, as the pressure decreases, the critical time occurs more slowly, and the risk of water inrush can be better controlled. However, the risk will be difficult to control and prone to water inrush hazards when the water pressure is larger.
4.3.2. Effect of Initial Effective Porosity
Since the formation process of each fault structure is different, the initial effective porosity of fault fillings is also different. By testing the seepage properties of fault fillings, it is possible to predict and prevent water inrush in advance.
Figure 14 plots the effect of different initial effective porosities on water flow rate. The results show that the effect of initial effective porosity on flow rate is obvious, and it is mainly concentrated in the second phase and third phase. The difference in flow rate is small in the initial phase under different initial effective porosity conditions. In the second phase, the difference in the effects of different porosity has become more significant. The main reason for this is that as the porosity increases, the coupling effect increases and the seepage characteristics evolve more strongly. By the final phase, the difference was significant and stable. However, the effect on the water inrush rate due to the difference in initial effective porosity is not easily noticeable in the initial phase. And it is easy to develop into an uncontrollable water inrush hazard in the fault with high initial effective porosity.
4.3.3. Effect of Initial Permeability
The initial permeability is also an important factor affecting the water flow rate. The initial permeability of the fault fillings can also be obtained by experiments to predict and prevent the water inrush hazards. Numerical simulations were carried out for different initial permeabilities under the same conditions, and the effect of different initial permeability on flow rate is shown in Figure 15.
As shown in Figure 15, for varying initial permeabilities, the evolution of both the water flow rate change and the water inrush critical time shows an obvious gap. As the initial permeability increases, the flow rate evolves more rapidly and intensely, and the critical time of water inrush occurs earlier. The difference in the effect of different initial permeability is significant throughout the process. For larger initial permeability, the evolution of water inrush is accelerated more intensely. Moreover, the larger initial permeability has a higher flow rate upper limit during the steady flow phase.
Initial water pressure, initial effective porosity, and initial permeability are all initial conditions for water inrush. When the values of these initial conditions are larger, the coupling effect tends to be significantly enhanced at the initial phase. At the same time, the evolution of seepage characteristics occurs earlier and more quickly. The result is that the water inrush time is greatly advanced and the upper limit of the flow rate is dramatically increased. However, aquifer water pressure, initial effective porosity, and initial permeability of fault fillings can all be measured in advance. If these conditions are combined with the above analysis and timely measures, it can effectively predict the critical time of water inrush and reduce the risk of water inrush hazards.
5. Conclusions
(i)In this paper, the coupling effect of seepage and particle migration in the process of water inrush is comprehensively considered, and a multifield coupling model for the evolution of characteristics of fault water inrush is established. The COMSOL software was used to numerically simulate the coupling model, and the simulation results of water inflow change were basically consistent with the field data, indicating the correctness of the evolution model of fault water inrush established in this paper(ii)According to the evolution process of seepage characteristics such as effective porosity, permeability, seepage velocity, and particle concentration, the evolution of water inrush process can be divided into three phases, including (i) low velocity seepage, (ii) drastic changes with substantial particle migration, and (iii) steady-state water flow. The seepage coupling will gradually increase with the evolution of seepage characteristics from the beginning. As the flow rate increases, the erosion of the water flow accelerates particle migration. This makes it easier for areas of large porosity and permeability to expand and connect to form cracks. Then, the overall permeability of the fault is enhanced, which in turn accelerates the evolution of the seepage velocity. This mutually reinforcing process continues to circulate, resulting in a rapid increase in the flow rate of water inrush. By the time the effective porosity has evolved to a maximum, the formation of the seepage channel and the evolution of the seepage characteristics have been stably completed(iii)Water flow rate increases with increasing initial water pressure, initial effective porosity, and initial permeability. When the values of these initial conditions are larger, the coupling effect tends to be significantly enhanced at the initial stage. At the same time, the evolution of seepage characteristics occurs earlier and more quickly. The result is that the water inrush time is greatly advanced and the upper limit of the water flow rate is dramatically increased(iv)The evolution model of seepage characteristics established in this paper can effectively predict fault water inrush. The critical time of water inrush by conditions such as aquifer water pressure, initial effective porosity, and initial permeability of fault fillings is predicted, so that proper measures can be put in place timely. It is expected to achieve the goal of reducing the risk of water inrush hazards and ensuring safe production
Data Availability
The data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
The project is supported by the National Natural Science Foundation of China (51404146), the China Postdoctoral Science Foundation (2015M572067), the Natural Science Foundation of Shandong Province (ZR201807060138), the Postdoctoral Innovation Project of Shandong Province (152799), the Key R&D Projects in Shandong Province (2015GSF120016), the Qingdao Postdoctoral Applied Research Project (2015203), the Anhui Province Key Laboratory of Modern Mining Engineering Foundation (KLMME12101), the SDUST Research Fund (2018TDJH102), and the Young Teachers’ Growth Program of Shandong Province.