1 Introduction

Gas disasters significantly affect mine safety. Due to increased mining depth and intensity in China, many low-permeability gas mines have been transformed into high-permeability gas or gas outburst mines, adversely affecting safety (Al-Rubaie and Ben Mahmud 2020; Du et al. 2020; Liu et al. 2011; Wang et al. 2013). Coal bed methane exploitation can prevent coal and gas outbursts and provide clean energy. An increase in the mining depth has intensified the in-situ stress, decreased the permeability of coal seams, and increased the difficulty of coalbed methane exploitation (Huang and Liu 2018; Liu et al. 2019; Wang et al. 2018). Therefore, improving the permeability and preventing gas disasters are crucial in coalbed methane development (He et al. 2012; Huang et al. 2015; Li et al. 2016; Wang et al. 2017b; Zhao et al. 2020; Zhou et al. 2021). Research has shown that hydraulic fracturing increases coal seam permeability, which is conducive to the exploitation of coalbed methane and the prevention of coal and gas outburst (Cheng et al. 2018; Wang et al. 2014; Zhang 2014).

Hydraulic fracturing of coal seams refers to injecting high-pressure liquid into the coal seam through the borehole using a high-pressure pump. The pressure overcomes the in-situ stress, tensile strength, and cohesion of coal seams, forming a fracture network that improves gas flow. Wang et al. (2014) used PFC2D software to analyze the effect of the macro mechanical properties on crack initiation and size during hydraulic fracturing. The results showed that the main fracture extended into the direction of maximum principal stress. Zhou et al. (2017) used a two-dimensional particle flow simulation program and a smooth joint model to simulate hydraulic fracturing. The authors studied the interaction between hydraulic fractures and natural fractures and compared the results with indoor experimental results and analysis model results. Wang et al. (2018) used the finite element method to study the fracture propagation mechanism of hydraulic fracturing in a coal seam with discontinuous natural fractures. The results showed that the hydraulic fracture network was spindle-shaped with a multi-level branch structure at a high stress difference. Secondary fractures accounted for a large proportion, representing an important part of the hydraulic fracture network in the coal seam. Wang et al. (2017a) used the particle flow code (PFC) to simulate the evolution of the overburden fissure field; the research results provided valuable guidance for gas drainage in coal mine fields. Lyu et al. (2020) analyzed the influence of natural fractures in coal seams on hydraulic fracturing performances. It was observed that the propagation of hydraulic fractures depended on the direction of natural fractures and the principal stress. Yuan et al. (2012) established a mathematical model of hydraulic fracturing in a low-permeability coal seam and analyzed the influence of the water injection pressure and other factors on the characteristics of fracture expansion. The results demonstrated that the fracture length increased linearly, and the fracture width increased exponentially with an increase in injection pressure. The fracture width was larger than the fracture length in the later fracturing stage. Zhao et al. (2019) studied the influence of the coal rock type (bright coal, semi-bright coal, semi-dark coal, and dark coal) and perforation location on hydraulic fracture propagation. It was found that the fracture morphology and proppant distribution were better for bright coal than dark coal.

The literature indicates that many scholars have investigated the fracture propagation type, fracturing influence range, and the influence of existing fractures using experiments and simulations (Adachi et al. 2007; Li et al. 2014; Liang et al. 2017). Many studies were conducted on the influence of coal seam hydraulic fracturing on the gas drainage effect. Fan et al. (2019) proposed a fully coupled mathematical model of hydraulic stress damage considering gas–water two-phase flow to simulate hydraulic fracturing to enhance underground gas drainage. Zhang (2014) established a two-phase three-dimensional seepage hydraulic fracturing model for dual porous media and conducted simulations using well test data from a basin in Western China. The results showed that hydraulic fracturing promoted the desorption and diffusion of coalbed methane and significantly improved its production. Zhang et al. (2018) conducted a hydraulic fracturing study of the Nantong Mine in the southeast of the Sichuan Basin. Field investigations showed that hydraulic fracturing significantly increased the methane extraction rate; it was more than 10 times higher than that of conventional boreholes. The field research of Huang et al. (2012) showed that hydraulic fracturing improved the downhole gas permeability, and the gas drainage capacity increased by 15 times. The field application results of pulse hydraulic fracturing (Xu et al. 2017) showed that the proportion of micropores decreased by 7.7%, the proportion of mesopores increased by 23.1%, and the proportion of macropores increased by 2.9%, significantly improving the permeability of the coalbed methane reservoir. Li et al. (2015) observed that the gas desorption index \(k_{1}\) of the driving face fell below the critical value after hydraulic fracturing. The gas drainage volumes of the fracturing boreholes and pilot boreholes were 3.32 times and 3.07 times higher, respectively, than that of the normal boreholes.

Scholars have investigated coalbed methane development and the impact of hydraulic fracturing on coal seams. However, few studies analyzed the anti-reflection effect of hydraulic fracturing, and no suitable method exists to predict the antireflection effect of coal seams. In this paper, PFC discrete element software is used for hydraulic fracturing simulation. The connected fractures are extracted and imported into COMSOL Multiphysics numerical analysis software for gas drainage simulation to evaluate the permeability enhancement caused by hydraulic fracturing of the coal seam. The research results can be used to predict the anti-reflection effect of hydraulic fracturing in underground coal mines.

2 Fracturing simulation with PFC

2.1 Coupled model for hydraulic fracturing

The PFC5.0 software was used to establish a discrete element coupled model to simulate hydraulic fracturing. The fluid was stored in the pore grid, as shown in Fig. 1. Fluid exchange occurs in the adjacent pore grid due to the fluid pressure difference. The coupling of the fluid and solid is achieved by altering the contact force to change the channel pore pressure. The pressure is changed by modifying the mechanical characteristics of the study area. The pore pressure exerts a force on the particles inside the pores.

Fig. 1
figure 1

Watershed and flow path in particle binding model

The flow rate of the fluid exchange can be expressed by the Hagen–Poiseuille equation (Hubbert and Willis 1957; Shimizu et al. 2011):

$$q = ka^{3} \frac{\vartriangle p}{L}$$
(1)

where \(q\) is the flow rate, \(a\) is the opening of the fluid channel, which is related to the normal force of the two particles, \(k\) is the permeability coefficient, \(\vartriangle p\) is the pressure difference between the two pore basins, and \(L\) is the length of the fluid channel.

When the bond between the two particles is broken, the opening \(a\) is (Shimizu et al. 2011):

$$a = a_{0} + \lambda (d - R_{1} - R_{2} )$$
(2)

where \(d\) is the distance between two particles, \(R_{1}\) and \(R_{2}\) are the radii of two particles, and \(\lambda\) is the dimensionless multiplier.

At time \(\vartriangle t\), the change in pore fluid pressure due to fluid flow is (Hubbert and Willis 1957; Shimizu et al. 2011):

$$\vartriangle p = \frac{{K_{\text{f}} }}{{V_{\text{d}} }}\left( {\sum {q\vartriangle t - \vartriangle V_{\text{d}} } } \right)$$
(3)

where \(K_{\text{f}}\) is the compressive modulus of the fluid, \(V_{\text{d}}\) is the "domain", i.e., the pore volume, \(\vartriangle V_{\text{d}}\) is the pore volume change, and \(\sum\limits_{{}}^{{}} q\) is the total flow of the fluid.

2.2 Parameter calibration

Parameter calibration is performed to ensure that the macro-mechanical properties measured in the laboratory correspond to the meso-mechanical properties in the simulation experiments. The PFC contact model used for rock and soil mechanics analysis is the linear parallel bonding model (Pb model). The microscopic parameters of the Pb model correspond to the macro-mechanical parameters, and the microscopic parameters are calibrated by a uniaxial compression test and a uniaxial tensile test. As shown in Table 1, a model with a width of 1 m and a height of 2 m was adopted for the simulation. The interval of the particle radius between the maximum and minimum values is uniform. The minimum radius is 0.01 m, the ratio of the maximum radius to the minimum radius is 1.6, and the porosity is 0.06.

Table 1 Simulation basic parameters

Uniaxial compression and uniaxial tension numerical experiments are used to calibrate the parameters. The calibration result is shown in Fig. 2, and the meso-mechanical parameters after calibration are listed in Table 2. These data are used to perform the uniaxial compression and uniaxial tension simulations to obtain the macroscopic mechanical properties corresponding to the meso-mechanical parameters. Table 3 indicates that the measured and simulated macroscopic mechanical properties are very close.

Fig. 2
figure 2

Numerical experiments: a Uniaxial compression b Uniaxial tension (red line indicates crack)

Table 2 Calibrated micromechanical parameters
Table 3 Comparison of measured and simulated macro mechanical properties

2.3 Analysis of fracturing effect

The numerical model is a two-dimensional model, as shown in Fig. 3. The length and width of the model are 5 m, the minimum and maximum radii of the particles are 0.010 m and 0.016 m, respectively, and the interval between the maximum radius and the minimum radius is uniform. The number of elements in the model is 43,457. The initial stress is simulated by adjusting the wall speed so that the horizontal stress of the model reaches 5 MPa and the vertical stress is 8 MPa. The water injection hole is located in the middle of the model and has a radius of 0.1 m. The fluid is injected at a rated pressure of 15 MPa. The mechanical model of the coal seam is a linear Pb model, and the parameter values are the micro-mechanical parameters obtained from the calibration.

Fig. 3
figure 3

PFC2D model before injection

Due to the long-term triaxial stress of the overlying strata, the coal seam has a certain strength. According to the stress conditions around the fracturing hole and the classical fracturing theory, many scholars have found that the fracture pressure is related to the horizontal effective stress, tensile strength, and pore pressure (Li et al. 2015). Therefore, the fracture pressure of the coal seam can be calculated as follows:

$$P_{\text{ini}} = 3{\overline{\sigma }}_{3} - {\overline{\sigma }}_{1} + P_{\text{pore}} +\sigma _{\rm t}$$
(4)

where \(\upsigma _{\rm t}\) is the tensile strength of the coal seam (MPa). \({\overline{\upsigma }}_{1}\) and \({\overline{\upsigma }}_{3}\) are the horizontal maximum and minimum principal stress (MPa), respectively. \(P_{\text{pore}}\) is the pore pressure.

The pressure near the injection hole is measured in a circle with a radius of 0.12 m. The relationship between the injection pressure and the number of cracks over time is shown in Fig. 4. The peak pressure of 7.42 MPa is called the burst pressure. The rupture pressure of the model calculated by Eq. (4) is 9 MPa, which is 1.58 MPa higher than the simulated rupture pressure of 7.42 MPa. The difference can be attributed to the excessively large initial pressure of the model and the short pressure accumulation time. Since the liquid is continuously injected, it accumulates at the injection orifice and the previously formed fracture. The pressure at the crack tip rises over time until it reaches the fracture pressure of the coal seam. New cracks appear in the coal seam, causing a pressure drop at the injection orifice. Therefore, the injection pressure curve shows fluctuations as the fracture expands. It is observed in Fig. 4 that the rate of fracture propagation changes from fast to slow over time. Moreover, the crack propagation is discontinuous, and the pressure at the tip of the crack is not sufficient to damage the coal body initially. The crack continues to expand after the pressure has reached the failure pressure.

Fig. 4
figure 4

Injection pressure and fracture growth curve

The fracturing effect at different time steps is shown in Fig. 5. The fractures expand in the horizontal and vertical directions but dominantly in the maximum principal stress direction. The propagation distance of the fracture is 2.38 m at 100,000 steps. As shown in Fig. 6a, the pressure of the fracturing fluid continues to increase at the fracture tip. Under the combined action of the fluid’s pressure and ground stress, the horizontal force of the particles at the fracture tip (the direction of the minimum principal stress) increases. As a result, the speed of horizontal crack propagation increases. Eventually, the particles move in the direction of the minimum principal stress, and the cracks continue to expand in this direction. As shown in Fig. 6b, the particles have all moved in the direction of the minimum principal stress, and the closer they are to the injection hole, the greater the displacement is, the more pronounced the crack propagation is.

Fig. 5
figure 5

Fracturing effect diagram at different time steps

Fig. 6
figure 6

Particle vector diagram under 100,000 fracturing steps: a velocity diagram b displacement diagram

3 Model establishment and simulation scheme

3.1 Model establishment

Based on the theory of poroelasticity and seepage mechanics (Xie et al. 2014), Zhang et al. (2019) established a fully coupled mathematical model of coal deformation and gas flow considering the deformation attributed to adsorption and desorption.

3.1.1 Control equation of coal seam deformation

According to the total strain equation of coal, the Langmuir volume strain equation \(\upvarepsilon _{\text{s}} =\upvarepsilon _{L} \frac{p}{{p + p_{L} }}\), the balance equation \(\upsigma _{ij,j} + f_{i} = 0\), and the Cauchy equation \(\upsigma _{ij,j} + f_{i} = 0\), the governing equation for coal deformation is obtained as:

$$\begin{aligned}&Gu_{i,kk} + \frac{G}{1 - 2\nu }u_{k,ki} - \upalpha p_{,i}\\ &\quad - K\upvarepsilon _{L} \frac{{p_{L} }}{{(p + p_{L} )^{2} }}p_{,i} + f_{i} = 0 \end{aligned}$$
(5)

where \(\upvarepsilon _{{{{ij}}}}\) is the component of the total strain tensor, \(G\) is the shear modulus of coal, \(K\) is the bulk modulus of coal, \(E\) is Young's modulus of coal, \(\nu\) is the Poisson's ratio of coal, and \(p\) is the gas pressure in the pores. The effective stress component is defined as \(\sigma _{{_{ij} }}^{\prime } =\sigma _{ij} +\sigma {{p}}\delta _{ij}\), \(\alpha = 1 - {{K}}/{{K}}_{{{s}}}\) is the Biot coefficient, \({{K}}_{{\text{s}}}\) is the bulk modulus of the coal particles, \(\delta _{ij}\) is the Kronecker number, \(u_{i}\) is the component of displacement, \(\upsigma _{ij}\) is the component of the total stress tensor, \(f_{i}\) is the component of the force on the object, \(\upvarepsilon _{L}\) is the Langmuir pressure constant, and \(p_{L}\) is the pore pressure. At this time, the measured volumetric strain is equal to \(0.5\upvarepsilon _{L}\) (Zhang et al. 2008).

3.1.2 Gas flow control equation

According to the ideal gas law, mass balance equation, and Darcy's law, the gas flow control equation is:

$$\left[ {\phi + \frac{{\uprho _{c} p_{\text{a}} V_{L} p_{L} }}{{(p + p_{L} )^{2} }}} \right]\frac{\partial p}{{\partial t}} + p\frac{{\partial\phi }}{\partial t} - \nabla .\left( {\frac{k}{\upmu }p\nabla p} \right) = Q_{s}$$
(6)

where \(p_{\text{a}}\) is then atmospheric pressure (101.325 kPa), \(k\) is the permeability, \(\mu\) is the dynamic viscosity of the gas, \(\upphi\) is the porosity, \(\rho _{\text{ga}}\) is the gas density under standard conditions, \(\rho _{\text{c}}\) is the coal density, \(\rho _{\text{g}}\) is the gas density, \(\nu _{\text{g}}\) is the Darcy velocity vector, \(Q_{\text{S}}\) is the gas source, \(t\) is the time, \(s\) and \(m\) are the gas contents of the free state gas and adsorbed gas, respectively (Mi et al. 2018).

3.1.3 Porosity and permeability models

It is assumed that the adsorption strain of coal is the same as the adsorption strain of the pore space. When the initial pressure is \(p_{0}\), the initial porosity is \(\upphi _{0}\), and the initial volumetric strain is zero. According to the site conditions, it is assumed that \(\upvarepsilon _{33}\) is the uniaxial strain direction, and the overburden load direction is unchanged. \(\upvarepsilon _{11}\) and \(\upvarepsilon _{22}\) are lateral strains and both are zero (Kumar et al. 2016; Mi et al. 2018; Pan and Connell 2012; Zhang et al. 2008). The permeability expression for the given uniaxial strain and a constant overburden load is:

$$\begin{aligned} k &= k_{0} \left\{ 1 + \frac{\alpha }{{\phi _{0} }}\left[ \frac{\alpha }{M}(p - p_{0} ) + \frac{{p - p_{0} }}{{K_{\text{S}} }}\right.\right.\\ &\quad\left.\left.+ \left( {\frac{K}{M} - 1} \right)\upvarepsilon _{L} p_{L} \frac{{(p - p_{0} )}}{{(p_{0} + p_{L} )(p + p_{L} )}} \right] \right\}^{3} \end{aligned}$$
(7)

where M is the limiting axial elastic modulus, \(M = E(1 -\nu )/(1 +\nu )(1 - 2\nu )\).

The effect of grain compression is considered for the rebound pressure of the fully coupled pore model. If the uniaxial strain and the overburden load remain unchanged, the expression of the rebound pressure is (Mi et al. 2018):

$$p_{\text{c}} = \left[ {\frac{{K_{\text{S}} (M - K)\upvarepsilon _{L} p_{L} }}{{\upalpha K_{\text{S}} + M}}} \right]^{1/2} - p_{L}$$
(8)

The partial derivative of the porosity \(\upphi\) with respect to time \(t\) is substituted into Eq. (6) to obtain the fully coupled control equation of gas flow under the influence of coal seam deformation:

$$\begin{aligned} &\left[ {\phi + \frac{{\rho _{\text{c}} p_{\text{a}} V_{L} p_{L} }}{{(p + p_{L} )^{2} }} + \frac{{(\alpha -\phi )p}}{{K_{\text{S}} }} - \frac{{(\alpha -\phi )\upvarepsilon _{L} p_{L} p}}{{(p + p_{L} )^{2} }}} \right]\frac{\partial p}{{\partial t}}\\ & \quad - \nabla \cdot \left( {\frac{k}{\mu }p\nabla p} \right) = Q_{\text{S}} - (\alpha -\phi )p\frac{{\partial\upvarepsilon _{v} }}{\partial t} \end{aligned}$$
(9)

Therefore, Eqs. (5), (6), (7), (8), and (9) define the fully coupled model of coal seam deformation and gas flow. The coupling relationship between coal seam deformation and gas flow is shown in Fig. 7.

Fig. 7
figure 7

Coupling relationship between coal seam deformation and gas flow

3.2 Simulation scheme and parameters

The above-mentioned control equations constitute a mathematical model describing coal deformation and gas flow. COMSOL Multiphysics software is used to solve the numerical model using the finite element method. Fractures are important flow channels during gas drainage; thus, different fracturing times result in different fracture network conditions and gas drainage effects. Different initial pore pressures also affect the gas drainage effects. Therefore, in this work, we analyze the effects of different fracturing times and different initial pore pressures on the gas drainage effect. The specific simulation schemes are shown in Table 4.

Table 4 Experimental conditions of different samples

First, the fracture propagation diagram for different fracturing time steps is imported into a CAD program, and the connected fractures are extracted, as shown in Fig. 8. The results are imported into COMSOL Multiphysics software. The calculation grid is densified to ensure the integrity of the fracture boundary. The final physical model of gas drainage is shown in Fig. 9. The size of the model is all 5 m × 5 m, the suction pressure is one atmosphere, and the boundary of the model is rolling support. The main parameters used in the numerical simulation are listed in Table 5.

Fig. 8
figure 8

CAD extraction of fracture propagation effect diagram of different fracturing time steps

Fig. 9
figure 9

Gas drainage physical model under different fracturing time steps

Table 5 Main parameters used in numerical simulation

4 Results and discussion

4.1 Fracture and non-fracture models

The established mathematical model is numerically solved using the parameters of the physical model. Figure 10a and b show the gas pressure distribution of the non-fracture model (S1) and the fracture model (S5) at different times of drainage. The results show that the gas pressure decreased rapidly in the first few days, but the decrease rate decreased slowly with an increase in the drainage time. The gas pressure of the fracture model decreased faster than that of the non-fracture model. The former shows differences in the distribution of the gas pressure in the horizontal and vertical directions. The gas pressure was lower in the vertical direction (main direction) of fracture propagation than in the horizontal direction. For the quantitative analysis, we used a detection point (1.5, 0) in the non-fracture model (because the model is symmetric in the X and Y directions, only one monitoring point is required) and two in the fracture model, point 1 (1.5, 0) and point 2 (0, 1.5), as shown in Fig. 11. The results show that the gas pressure trend of the two models was the same. The gas pressure dropped quickly at the beginning of the drainage, and the gas pressure decline rate rapidly decayed and approached zero as the drainage time increased. The gas pressure of the fracture model decreased faster than that of the non-fracture model. In the fracture model, the gas pressure decreased faster at monitoring point 2 than at monitoring point 1, because monitoring point 2 was located in the main direction of fracture propagation. After 50 d of drainage, the gas pressure at the monitoring point of the non-fracture model dropped from 2 to 0.48 MPa, that at monitoring point 1 of the fracture model dropped to 0.22 MPa, and that at monitoring point 2 dropped to 0.18 MPa. According to national standards, the gas pressure should be reduced to below 0.74 MPa to eliminate the gas outburst risk. It required 25 d for the gas pressure at the monitoring point of the non-fracture model to fall below the standard. In contrast, it took only 8 d and 6 d for the fracture model at monitoring point 1 and monitoring point 2, respectively. The gas pressure curve obtained by Jia et al. (2020) at different distances from the fracturing hole is shown in Fig. 12. The results show that the gas pressure decreased rapidly in the initial stage of drainage and then gradually approached 0 with an increase in the drainage time. The further the distance from the fracturing hole, the slower the gas pressure drop was. This result shows that the anti-reflection effect of hydraulic fracturing becomes worse with the increasing distance from the fracturing hole, which is consistent with the simulation results in this paper.

Fig. 10
figure 10

The distribution of coal seam gas pressure at different times: a Unfractured model b Fractured model

Fig. 11
figure 11

Gas pressure variation curve with time at different monitoring points

Fig. 12
figure 12

Variation curve of gas pressure at different intervals from fracturing hole

Figure 13 shows the permeability increase rate \((k - k_{0} )/k_{0}\) of the non-fracture model and the fracture model at different times of drainage. The results show that as the extraction time increases, the permeability increase rate is higher for the fracture model than the non-fracture model. However, the increase rate of permeability in the fracture model is different in the horizontal and the vertical directions, and the increase rate of the permeability is slightly higher in the main direction of fracture propagation. A detection line was created between points (0, 0) and (2.5, 0) of the non-fracture model and the fracture model to quantify the rate of change of permeability over time, as shown in Figs. 14 and 15. The results show that the further from the center point, the smaller the increase rate of permeability was. The decrease rate gradually approached 0, changing from fast to slow, and the permeability increase rate was higher for the fracture model than the non-fracture model. As the drainage time increased, the increase rate of permeability rose significantly. After 20 d of drainage, the permeability at the center and edge of the non-fracture (fracture) model increased by 2.1% and 1.1% (2.27% and 1.7%), respectively.

Fig. 13
figure 13

The distribution of coal seam permeability increase rate at different times: a Unfractured model b Fractured model

Fig. 14
figure 14

Permeability increase rate change curve of the unfractured model

Fig. 15
figure 15

Permeability increase rate change curve of fracturing model

4.2 Different fracturing times

Figure 16 shows the increase rate of gas pressure and permeability at different time steps (S2–S5) of fracturing for different gas drainage times. As the drainage time increased, the gas pressure declined, and the permeability increase rate increased. The gas pressure decrease rate and permeability increase rate were higher in the direction of maximum principal stress (main direction of fracture propagation) than the direction of minimum principal stress for all fracturing time steps. As the fracturing time increased, this difference also increased. Two detection points were used in the S2–S5 model: monitoring point 1 (1.5, 0) and detection point 2 (0, 1.5). The gas pressure at different fracturing time steps is shown in Fig. 17. The trends of the four models were the same. The gas pressure dropped rapidly in the early stage of drainage. As the drainage time increased, the gas pressure decline rate decreased sharply and tended to zero. The longer the fracturing time, the more complex the fracture network was, and the faster the gas pressure dropped. Moreover, as the fracturing time increased, the difference between the horizontal and vertical gas pressure decrease rates increased, and the rate was significantly higher in the direction of the maximum principal stress than the direction of the minimum principal stress. It took 13 d for the gas pressure to drop below 0.74 MPa for model S2 at monitoring points 1 and 2, 11 d and 10 d for model S3, 10 d and 9 d for model S4, and 8 d and 9 d for model S5. Models S3 and S5 required 11 d and 8 d, respectively, and the fracturing time of S5 was more than twice that of S3. The result shows that the influence of fracturing on the gas drainage effect decreases over time. The permeability increase rate at different fracturing time steps is shown in Fig. 18. As the drainage time increased, the rate of increase in permeability changed from fast to slow. Since test point 2 was located in the main fracture propagation direction, the increase in permeability of test point 2 was always greater than that of test point 1. The difference because more pronounced as the fracturing time increased. However, the difference decreased as the extraction time increased. The fracture radius and fracture width near the fracturing orifice during the fracturing process measured by Bao et al. (2021) are shown in Fig. 19. The fracture radius and width increased rapidly at the initial fracturing stage, but the rate of increase decreased sharply over time. This result shows that the influence of fracturing on the gas drainage effect decreased over time, which is consistent with the simulation results in this paper.

Fig. 16
figure 16

The distribution of different fracturing time steps after 20 days of pumping: a Gas pressure b Permeability increase rate

Fig. 17
figure 17

Gas pressure variation curve at different fracturing time steps

Fig. 18
figure 18

Variation curve of permeability increase rate at different fracturing time steps

Fig. 19
figure 19

Variation curve of fracture radius and fracture width with fracturing time

4.3 Different initial pore pressures

The same detection points were used in the S5–S8 models. Figure 20 shows the gas pressure of the models with different initial pore pressures (S5–S8). The results show that the gas pressure decline rate changes from fast to slow, and the gas pressure difference between point 1 and point 2 decreased over time. The higher the initial pore pressure, the faster the gas pressure dropped. It approached the orifice pressure as the extraction time increased. It required 8 d (point 1) and 6 d (point 2) for the gas pressure to drop below 0.74 MPa for model S5, 10 d and 7 d for model S6, 10 d and 8 d for model S7, and 11 d and 8 d for model S8. The initial pore pressure of model S8 was 1.5 times higher than that of model S5. However, after fracturing, the gas pressure dropped below the standard due to the drainage effect, and it only took 3 more days for extraction. Figure 21 shows the permeability increase rate of the models with different initial pore pressures. The trends of the different models were consistent. The permeability increased rapidly at the beginning of the drainage, but the increase rate decayed rapidly and approached zero as the drainage time increased. The greater the initial pore pressure, the greater the increase in permeability was. However, the influence of the initial pore pressure on the increase in permeability gradually decreased. For a quantitative characterization of the decline rate of the coal seam gas pressure, the gas pressure decline rate \(D_{\text{p}}\) was defined by the ratio of the gas pressure decline (difference between adsorbed gas pressure \(p_{0}\) and residual gas pressure \(p_{t}\)) and the adsorbed gas pressure at the time of gas extraction \(t\):

$$D_{\text{p}} = \frac{{p_{0} - p_{t} }}{{p_{0} }} \times 100\%$$
(10)
Fig. 20
figure 20

Gas pressure variation curves of different initial pore pressure models

Fig. 21
figure 21

Permeability increase rate change curve of different initial pore pressure models

Zhang et al. (2014) obtained the gas pressure and gas pressure decline rate over time from models with different initial pore pressures, as shown in Fig. 22. The experimental results show that the gas pressure decreased from fast to slow as the drainage time increased. The higher the initial pore pressure, the faster the gas pressure decreased, which is consistent with the simulation results in this paper.

Fig. 22
figure 22

Variation curves of gas pressure and gas pressure drop rate of different initial pore pressure models

4.4 Application prospects

This study analyzed the influence of hydraulic fracturing on coal seam permeability for different fracturing times and different initial pore pressures. The research results showed that hydraulic fracturing could increase the permeability of the coal seam and improve the gas drainage effect. However, as the fracturing time increased, the impact of fracturing on the gas drainage effect declined. The research results can be used to predict the increase in coal seam permeability due to hydraulic fracturing and adjust the fracturing time accordingly. The higher the initial gas pressure, the faster the gas pressure dropped, and the greater the increase in permeability was. The results can be used to guide hydraulic fracturing to enhance permeability in mines with different gas contents.

Due to in-situ stress, the main direction of fracture propagation is in the direction of maximum principal stress. Therefore, a difference exists in the permeability change between the horizontal and vertical directions, and this difference increases with the fracturing time. We will focus on eliminating or exploiting this difference in our next study.

5 Conclusions

In this work, a coupled mathematical model of the interaction between coal seam deformation and gas flow was established, and simulations were conducted to determine the effects of hydraulic fracturing on gas drainage. The research conclusions can be summarized as follows:

  1. (1)

    The fracture growth rate changed from fast to slow as the fracturing time increased and crack propagation was discontinuous. When the pressure at the crack tip was insufficient to damage the coal, no new cracks occurred, and when the pressure was equal to the breaking pressure, the cracks continued to expand. The particles moved in the direction of the minimum principal stress, and the main fractures expanded in the direction of the maximum principal stress. The closer to the injection hole, the greater the displacement was, and the more adequate the fracture expansion was.

  2. (2)

    The gas pressure trends were the same for the fracture model and the non-fracture model. As the drainage time increased, the decline rate of gas pressure changed from fast to slow. The gas pressure decrease rate and permeability increase rate of the fracture model were significantly higher than those of the non-fracture model. However, the change in the gas pressure and permeability was different in the horizontal and vertical directions.

  3. (3)

    The longer the fracturing time, the more complex the fracture network was, the faster the gas pressure decreased, and the permeability increased, but the greater the difference was between the horizontal and vertical directions. However, this difference gradually decreased with an increase in the extraction time. The fracturing effect on the gas drainage decreased over time.

  4. (4)

    The greater the initial pore pressure, the faster the gas pressure dropped. It approached the orifice pressure as the extraction time increased. The increasing trend of the permeability was the same for the S5–S8 models. The increase rate of permeability decreased exponentially with the increase in the extraction time. The greater the initial pore pressure, the larger the permeability increase was. However, the influence of the initial pore pressure on the rise in permeability gradually declined over time.