Gas seepage laws based on dual porosity and dual permeability: Numerical simulation and coalbed methane extraction practice

Coalbed methane (CBM) is the very important unconventional energy resource. Enhancing extraction is the main method to improve the utilization of CBM and prevent coal mine gas disasters. To reveal the laws of gas migration during the underground gas extraction, the paper established the solid‐gas coupling model, regarding coal as the homogeneous elastic medium with dual pore‐fracture structure and dual permeability and considering the gas dynamic diffusion coefficient. Then, the regulations of gas migration in the in‐seam borehole were simulated by the COMSOL Multiphysics software and verified via the field trials of gas extraction. The results revealed that the permeability increased with time due to the coupling of matrix shrinkage and volume compression of coal, in which the matrix shrinkage played a leading role. The gas seepage velocity at the observation points can be divided into three stages: the rapid rise stage, the stable decline stage, and the stable and invariant stage. The gas extraction rate continued decreasing with time and finally tended to a fixed value. In the coal seam #4 of Zhongxing coal mine, the measured flow rates of three in‐seam drilling holes were in good agreement with the simulation results, which verified the correctness of the theoretical coupling model.

lations of gas migration in the in-seam borehole were simulated by the COMSOL Multiphysics software and verified via the field trials of gas extraction. The results revealed that the permeability increased with time due to the coupling of matrix shrinkage and volume compression of coal, in which the matrix shrinkage played a leading role. The gas seepage velocity at the observation points can be divided into three stages: the rapid rise stage, the stable decline stage, and the stable and invariant stage. The gas extraction rate continued decreasing with time and finally tended to a fixed value. In the coal seam #4 of Zhongxing coal mine, the measured flow rates of three in-seam drilling holes were in good agreement with the simulation results, which verified the correctness of the theoretical coupling model.

K E Y W O R D S
coalbed methane, dual porosity/dual permeability, dynamic diffusion, numerical simulation, solidgas coupling model extraction of coal seam can not only guarantee the safe production of high-gas mines, but also improve the use of CBM, and it has significant safety, and economic and environmental benefits to reduce greenhouse gas emissions. 14,15 Coal is a complex porous medium, and coal permeability is regarded as a key factor that characterizes the ability of fluids to flow through coal reservoirs. 5,16,17 Establishing a reasonable gas seepage model and revealing the regulations of gas migration in coal seam can provide a reasonable theoretical basis for the prediction of gas extraction quantity and design of gas extraction engineering. Many researchers have done a lot of research on coal seam gas seepage model. Zhao et al 18 considered that the change in ground stress and gas pressure would have an impact on the effective stress of coal, established a nonlinear coupled mathematical model of solid deformation and gas seepage without considering gas desorption expansion, and carried out numerical simulation of gas extraction in coal seams. Connell et al 19 established a gas extraction coupling model for the influence of effective stress, pore pressure, and adsorption/ desorption effect on coal seam permeability. Zhou and Wang et al 20 established a gas-solid coupling theoretical model considering the impact of geo-stress, gas pressure, and mining depth, and carried out simulation research by using COMSOL Multiphysics software. Based on the multifield coupling theory, Zhu et al 21 established a gas-solid coupling model of coal and gas considering the Klinkenberg effect. Valliappan et al 22 established a coupling model considering the influence of adsorbed methane diffusion from the matrix to pore based on the dual-porosity/single-permeability model and described the relationship among particle deformation, gas flow, and gas pressure in the coal seam. Liu et al 23 revised the P-M permeability model and established a gas-solid coupling model of coal and gas considering the adsorption/desorption effect, Klinkenberg effect, and coal skeleton compression effect. Liu et al 24 based on the P-M permeability model introduced the dynamic diffusion coefficient of gas, established the gas-solid coupling model of coal and gas with dual porosity/single permeability, and numerically solved the model by using COMSOL Multiphysics software to study the spatial and temporal evolution laws of permeability under different diffusion coefficients and the distribution characteristics of gas pressure. Coal is a kind of actually having the dual pore-fracture structure and has dual-permeability homogeneous elastic medium. But at present, a lot of scholars of coal seam gas seepage model have established in a single porosity/ single permeability or dual porosity/single permeability under the premise conditions, and do not conform to the dual-porosity/dual-permeability characteristics, thus unable to fully understand the laws of gas migration in the process of extraction.
On the basis of previous studies, this paper established the solid-gas coupling model, regarding coal as the homogeneous elastic medium with dual pore-fracture structure and dual permeability, combining the Klinkenberg effect and matrix shrinkage effect, and introducing the dynamic gas diffusion coefficient. Then, the laws of gas migration were simulated and analyzed, and verified via the field of gas extraction practice, in order to provide a reasonable theoretical basis for the gas extraction design.

| Fundamental assumptions
1. The coal body is regarded as a homogeneous isotropic elastic medium with dual pore-fracture structure and dual permeability, 25 as shown in Figure 1. 24 2. The free gas in the coal is an ideal gas, which occurs in adsorptive and free states and is transported in the matrix pores and fractures. 3. The influence of gas body force and coal seam temperature on gas extraction are not considered. 4. Gas flows into the fractures from the matrix in two ways, which are seepage and diffusion. The diffusion coefficient is not constant, and the seepage follows Darcy's law. 26,27 5. The coal is in the stage of linear elastic deformation, obeyed the generalized Hooke's law, and conformed to the small deformation hypothesis.

| Evolution model of porosity and permeability of coal
Based on the above assumptions, the coal matrix porosity and fracture porosity models were constructed, respectively, and F I G U R E 1 Physical model for coal containing methane and methane migration processes then, the coal matrix permeability and fracture permeability evolution models were established. The sorption-induced volumetric strain fitted onto the Langmuir-type curves and has been verified through experiments [28][29][30] and can be calculated by the Langmuir-type equation.
Matrix porosity can be defined as follows 31 : where where b is porosity of coal matrix; b0 is the initial porosity of coal matrix; b is the pore Biot coefficient is the pore stress variable of coal matrix; Q 0 is the initial porosity strain of coal matrix; V is the volumetric strain; V0 is the initial volumetric strain; K is the volume modulus of coal, MP a ; K s is the bulk modulus of coal skeleton, MP a ; L is a constant representing the volumetric strain at infinite pore pressure and the Langmuir pressure constant; P L is the Langmuir pressure constant and represents the pore pressure at which the measured volumetric strain is equal to 0.5 L , MP a ; P b is the gas pressure of coal matrix, MP a ; and P b0 is the initial pressure of coal matrix gas, MP a .
The permeability and porosity of coal matrix satisfy the cubic law, 32 which is defined as follows: where k b is the permeability of coal matrix and k b0 is the initial permeability of coal matrix.
By substituting Equation (1) into Equation (4), it can be obtained that: Fracture porosity can be defined as 33 : where b is the initial fracture aperture, m; Δb is the variation in aperture of fracture; L f is modified fracture stiffness, L f = aL n , MP a ; L n is the stiffness of fracture,MP a ; a is the matrix size, m; and ΔP b is the variation of coal matrix gas pressure, MP a .
It is known that there is a cubic law of permeability and porosity of coal fracture, 32 and the following conclusions are obtained: Klinkenberg effect: When gas flows in porous media, it has compressible and effective permeability is a function of gas pressure. When the average free path of gas molecules is equal to the pore size of porous media, the velocity of gas molecules near the surface of pore wall is not 0. 34 Because coal is a dual medium structure with pores and fractures, and the permeability of deep coal is low, the flow of gas in coal will have a significant Klinkenberg effect.
Considering the Klinkenberg effect comprehensively, the fracture permeability can be obtained as follows: where c is the Klinkenberg factor and P f is the gas pressure of coal fracture, MP a .

| Seepage and diffusion equations of gas in pores
The existing forms of gas in coal matrix include adsorption gas and free gas; therefore, the mass m of coal matrix gas per unit volume can be calculated by the Langmuir-type Equation 35 : where b is the gas density of coal matrix, kg/m 3 ; g is the gas density under standard condition, kg/m 3 ; c is the coal skeleton density, kg/m 3 ; and V L is the Langmuir volume constant and represents the volume of adsorbed gas per unit mass of coal saturated with adsorbed gas at a given temperature, m 3 /kg.
The diffusion coefficient is an important parameter reflecting the diffusion resistance of the medium. The paper 27 established a dynamic diffusion coefficient model with high accuracy and negative exponential change with time based on the analysis of the gas diffusion mechanism and the use of data fitting methods. The dynamic diffusion coefficient d t is as follows: where d t is the dynamic diffusion coefficient, m 2 ∕s; d 0 is the initial gas diffusion coefficient, m 2 ∕s; is the attenuation coefficient, s − 1 ; and t is the time, s. Free gas in coal matrix pore enters the fracture through diffusion and seepage. The diffusion coefficient is not constant, and seepage follows Darcy's law. According to the laws of gas migration, the law of conservation of mass, and the introduction of the dynamic gas diffusion coefficient, the seepage and diffusion field equations of gas in the matrix pores can be obtained 24 : where M g is the molar mass of methane, kg/m 3 ; R is the universal gas constant, J/mol K; T is the coal seam temperature, K; is the dynamic viscosity of gas,Pa s; s is the density of coal matrix, kg/m 3 ; is the matrix shape factor, = 3 2 a 2 , m -2 ; P s is the standard atmospheric pressure, 101.325kP a ; T s is the standard temperature, 273.15 K; and ∇ is a Hamiltonian operator.

| Seepage field equation of fractured gas
There is a mass exchange between the gas in the matrix and the fracture. The gas in the matrix is equivalent to an internal quality source of the fracture. The variation of the gas in the fracture is equal to the difference between the gas flowing into the fracture and the gas flowing into the extraction borehole. It can be expressed as follows:

| Governing equation of stress field
In this paper, based on the modified effective stress principle of dual-porosity medium, combined with gas adsorption/desorption, according to the hypothesis (5), the modified coal deformation equation based on the Navier form is as follows 35 : where u i,ij is displacement components (i, j = 1, 2, 3); G is the shear modulus of coal containing gas, MP a ; f is the Biot coefficient of fracture, f = 1 − K aL n ; and F i is the body force, MP a .
According to the permeability evolution models, the change in coal seam gas pressure (matrix gas pressure P b and fracture gas pressure P f ) was related to permeability. According to the cross-coupling model, the dynamic changes in matrix permeability k b and fracture permeability k f depend on the dynamic changes in coal matrix porosity b and fracture porosity f , respectively. With the progress of gas extraction, due to the decrease in gas pressure, and the increase in effective stress, and the compression of coal skeleton, and the shrinkage of coal matrix, leading to the corresponding changes of porosity b and fracture porosity f . In addition, the change in fracture permeability f is not only affected by the fracture porosity b , but also affected by the change in fracture gas pressure P f .
Gas migration in coal seam includes diffusion and seepage, and there is mass exchange between the two migration fields. The matrix gas is involved in diffusion and seepage after desorption, and the seepage results in the fracture have a reverse impact on the diffusion and mass transfer rates of the coal matrix after being affected, and they are boundary conditions for each other. Therefore, the mass exchange from the coal matrix to the fracture is called positive mass source, while the mass exchange from the fracture to the coal matrix is called negative mass source. According to the effective stress principle of dual media, the gas migration field will have a certain impact on the deformation field of the coal body. When the gas pressure balance in coal body is broken, the pressure gradient can be equivalent to the body force.

| Numerical model and parameters
The numerical simulation in this paper takes the coal seam #4 of Zhongxing coal mine in Fenxi as the actual engineering background to reveal the temporal and spatial evolution characteristics and changing laws of various parameters in the process of extraction. The source of the parameters used in the numerical simulation is the engineering measurement results, and it refers to other recently published papers. 24,35,36 The simulation parameters used in this paper are shown in Table 1.
In this paper, a single borehole in-seam is taken as the gas extraction model. COMSOL Multiphysics software was used to construct a 2D simplified horizontal section geometric model of the coal seam, and the model meshing adopts the built-in free triangle mesh of the software. The geometric model parameters are as follows: coal seam height W = 5 m and length D = 50 m; the diameter of the borehole is d = 80 mm, and the borehole is located in the center of the model, and observation points A 1 , A 2 , A 3 , and A 4 are set up at the horizontal distances of 0.5, 1, 1.5, and 2 m from the right side of the borehole. The solution time is 180 days.
As shown in Figure 3, coal is a typical dual-porosity system, in which two pressures are defined at each point: gas pressure P f in fractures and gas pressure P b in matrix pores. P b is defined as the "virtual" pressure in the matrix pores because the gas pressure in the matrix pores is difficult to be discussed clearly. Most of the analysis models assume that the system has reached the equilibrium state in the initial state, that is, the gas pressure in the matrix pore and the gas pressure in the fracture are equal to the gas pressure in the coal reservoir. 37 Therefore, the initial pore gas pressure and fracture gas pressure of the model are both 0.8 MPa. In addition, the relevant Chinese regulations

Parameter Value
Initial porosity for coal matrix, b0 0.05 Initial porosity for coal fracture, f0 0.001 Initial permeability of coal matrix, k b0 (mD) 1 × 10 −5 Initial fracture permeability of coal, k f0 (mD) 0.1 Initial fracture aperture, b(m) 35  stipulate that the negative pressure shall not be less than 13 kPa, so the fracture and pore pressure at the borehole boundary are set to be 0.088 MPa. For diffusion and seepage, referring to the paper, 21 the no-flow boundary condition is applied to the four sides of the model. The actual situation meets the uniaxial compression condition, and the top is the pressure boundary, which bears the overburden pressure of 10 MPa. In fact, all the boundaries of coal body should be set as stress boundary, which is most in line with the site conditions. However, it is found that if all the boundaries of the model are set as stress boundaries, the model cannot be calculated. Roller support means that there is no displacement in the normal direction, and the left and right sides of the model are set as roller support. 24 The bottom of the model is set to a fixed boundary. Figure 4 shows the distribution of matrix gas pressure nephogram at 60, 120, and 180 days after gas extraction. The variation of gas pressure can be used as a characterization factor of gas flow and an important evaluation parameter in the actual gas extraction process. Most of the gas in the coal seam is stored in coal matrix in the form of adsorption state, so the distribution of matrix pore pressure P b is discussed first. According to Figure 4, the cloud chart of matrix gas pressure distribution near the borehole is elliptical funnel-shaped.

| Variation law of matrix pore gas pressure
The gas pressure in the center area of the funnel is the lowest, and the gas pressure value at the outermost side of the funnel, that is, the farthest distance from the drilling center, is the largest. The radius of the funnel gradually expands with the increase in extraction time, and the gas pressure decreases with the increase in time. The maximum gas pressure of the whole coal seam matrix when the extraction reaches 180 days is only 0.526 MPa. Compared with the initial matrix pore gas pressure value of 0.8 MPa, the pressure value decreased significantly. The gas flow of dual-porosity/dual-permeability model is divided into two stages: The first stage is the process that gas in the matrix enters the fracture through diffusion and seepage. In this stage, the permeability is very small, the seepage resistance is particularly large, and the pressure propagation and decline are very slow; the second stage is that the gas in the fracture enters into the extraction borehole through seepage, with less seepage resistance and faster pressure drop. Figure 5 shows the distribution of matrix pore gas pressure and the distance from the borehole at different extraction time. As shown in Figure 5, the closer to the extraction borehole, the lower the gas pressure in matrix pore, and with the increase in the distance from the borehole, the change in gradient of gas pressure gradually decreases. The reason is that the closer to the borehole center, the more obvious the pressure relief effect of the local coal seam, the greater the coal seam permeability, and the more obvious the gas pressure drop. The gas pressure gradient shows that in the early stage of gas extraction, the gas extraction volume is relatively large.  Figure 6 shows the distribution of fracture permeability at different extraction time. As shown in Figure 6, the permeability near the borehole is much higher than the permeability farther away from the borehole, and the rate of decrease gradually slows down. The more obvious the effect of pressure relief near the borehole, the greater the permeability, so that the permeability of the area near the borehole is greater than that at a distance from the borehole. Figure 7 shows the variation of the fracture permeability with time at observation points A 1 (25.5,2.5), A 2 (26,2.5), A 3 (26.5,2.5), and A 4 (27,2.5). As shown in Figure 7, the fracture permeability at the observation points increases with the increase in time. Based on the solid-gas coupling model established above, it can be seen that gas extraction reduces the pressure of coal seam gas, and under the condition that the force of overlying strata is almost unchanged, the effective stress on coal body increases continuously, which leads to the compression of coal skeleton, porosity, and permeability; the gradual decrease in pore pressure leads to the shrinkage of the coal matrix and the increase in permeability. The coupling results of the coal matrix shrinkage and the coal skeleton compression lead to the dynamic change in permeability. The change in fracture permeability is closely related to the change in gas pressure. The lower the gas pressure, the more obvious the matrix shrinkage effect, which leads to the greater permeability. As the pressure continues decreasing, the degree of matrix shrinkage continues to increase. Under the influence of borehole extraction, although the coal skeleton is compressed to some extent under the action of stress field, the matrix shrinkage caused by desorption and diffusion of matrix gas is more obvious and far greater than that caused by stress field, so that the fracture permeability increases with the increase in time. At the initial stage, the fracture pressure decreases rapidly, the matrix gas flows rapidly under the action of large fracture pressure difference, the matrix adsorbed gas rapidly desorbs and flows, the pressure also drops rapidly, the matrix shrinkage effect is enhanced, and the permeability growth rate is also fast. Figure 8 shows the gas seepage velocity in the fractures of observation points A 1 , A 2 , A 3 , and A 4 and the gas extraction velocity in the borehole over time. It can be seen from Figure 8 that the seepage velocity at observation points rises rapidly at the initial stage and reaches the peak value, and subsequently starts to decline steadily and finally tend to be stable. The variation of seepage velocity with time can be roughly divided into three stages: the rapid rise stage, the stable decline stage, and the stable and invariant stage. At the initial stage, due to the large pressure difference between the fracture and the extraction borehole, the free gas in the fracture quickly enters the extraction borehole, which makes the seepage velocity of the fracture rise rapidly. With the extraction going on, the free gas content in the coal seam fracture decreases. Although there are matrix gas desorption/diffusion and seepage supplement, the supplement speed is less than that from the fracture to the extraction borehole. In the later stage of extraction, the difference between the matrix gas pressure and the fracture gas pressure tends to a fixed value, which makes the fracture seepage velocity maintain a certain value, thus entering the stable seepage stage. As can be seen from the gas extraction line of the borehole in Figure 8, at the beginning, the gas extraction rate decreased by a large gradient, and the rate of gas extraction decreased rapidly from 42.41 to 3.40 m 3 /d. Later, it decreased slowly with the growth of extraction time, and the rate of reduction was much smaller than that at the beginning. Finally, it gradually tended to be stable, and the stable value was about 1.7 m 3 /d.

EXTRACTION PRACTICE
The Zhongxing coal mine is located in the central part of Shanxi Province. The main coal seams are #2, #4, and #5. About 6m below coal seam #2, the coal seams #4 and #5 in the central and western regions are combined. The combined area is called the #4 coal seam, and most of the coal seams are mineable. The average thickness is 4.66 m, and the buried depth of the coal seam is 344.9-704.1 m. The original gas content of coal seam #4 is 5.43-7.05 m 3 /t, and the maximum gas pressure is 0.82 MPa, which is a difficult-to-drain coal seam. #4 and #5 coal seams slate are mudstones and carbonaceous mudstone, and floor lithology is sandy mudstone and mudstone. The single longwall mining method, retreating mining, and roof management by total caving method are adopted for the coal seams mined in the mine. The location of the experimental area is shown in Figure 9.
In order to verify the correctness of the theoretical coupling model, three in-seam drilling holes were constructed in the coal seam #4. The boreholes diameter was 80 mm, the length was 120 m, and the borehole spacing was 6 m, as shown in Figure 10. The attenuation laws of gas flow in in-seam boreholes were continuously measured on site, as shown in Figure 11.
It can be seen from gas extraction purity and mixing quantity of borehole monitored in Figure 11 that the gas extraction purity rates of the three in-seam boreholes decreased rapidly from 47.52, 39.36, and 55.68 m 3 /d at the beginning to 7.65, 4.86, and 5.13 m 3 /d, respectively, and then decreased gradually and finally reached the stable value of 1.44, 1.49, and 1.69 m 3 /d.
Through the comparison between the gas extraction pure quantity of the monitored boreholes in Figure 11 and the numerical simulation. It can be concluded that during the whole extraction time, the values of the gas extraction rates of the three monitored boreholes fluctuate up and down the simulated gas extraction rate. The simulation result curve has several intersections with the three boreholes gas extraction pure quantity curves. At the beginning, the difference between the monitored values and the simulated gas extraction rate was large, and then, it gradually decreased with the growth of time, and the last four curves were close to coincidence. It can be concluded that the gas extraction rates of the monitored boreholes are not too different from that of the numerical simulation results and have a good consistency, thus verifying the reliability of the theoretical coupling model.

| CONCLUSIONS
In this paper, we established the solid-gas coupling model with dual porosity/dual permeability in gas-bearing coal seams, and the laws of gas migration in the in-seam borehole under different extraction conditions were simulated and analyzed, and the field trials of the underground gas extraction were used to verify it. The major conclusions are summarized as follows: Established the solid-gas coupling model of dual porosity/dual permeability in gas-bearing coal seam, regarding coal as the homogeneous elastic medium with dual pore-fracture structure and dual permeability and considering gas adsorption/desorption, Klinkenberg effect, matrix shrinkage effect, and introducing the gas dynamic diffusion coefficient. The gas flow and mass transfer affect the matrix gas pressure P b and the fracture gas pressure P f , and the two pressures also affect the permeability and diffusion, thus affecting the entire dual-porosity/dual-permeability process. 2. Numerical simulation models were established to study the dynamic evolution laws of matrix gas pressure, fracture permeability, and gas drainage rate under different drainage time conditions. The evolution of permeability is the competition result of two opposite effects: the matrix shrinkage effect and the coal skeleton compression.
The gas seepage velocity at the observation point can be divided into three stages: the rapid rise stage, the stable decline stage, and the stable and invariant stage. The gas extraction rate of the coal seam is gradually decreased with time. The rate of decrease was relatively fast at the beginning, gradually slowed down with the increase in time, and finally tends to a stable value of 1.7 m 3 /d. 3. In order to verify the correctness of the model, three inseam drilling holes were constructed in the coal seam #4 to monitor the attenuation laws of gas flow in the boreholes. The difference between the initial field monitored results and the numerical simulation was high, and the difference gradually decreased with the increase in time and finally tended to coincide. This showed that the measured flow rates of three in-seam drilling holes were in good agreement with the simulation results, which verified the correctness of the theoretical coupling model.