Non-Darcy Flow in Hydraulic Flushing Hole Enlargement-Enhanced Gas Drainage: Does It Really Matter?

In recent years, the non-Darcy flow has seen a significant increase in interest in conventional and unconventional gas developments. The nonlinear behavior can be described by adding a quadratic term of the velocity with the permeabilitydependent β factor in Darcy’s law. The methods for improving permeability such as the hydraulic flushing hole enlargement (HFHE) technique can remarkably enhance the coal permeability thus changing the β factor. However, few studies have been done previously to seek how this non-Darcy flow impacts the permeability-enhanced gas drainage. In this study, a fully coupled coal deformation and damage, the non-Darcy flow of free methane gas in the fractures, and Fickian diffusion of adsorbed methane in the coal matrix model were developed. The Klinkenberg effect has been taken into account as a reference object. A series of simulation scenarios were carried out to evaluate and compare the influences of the non-Darcy effect, gas-coal interaction effect, and Klinkenberg effect on the HFHE-enhanced gas drainage. Results of the simulation illustrate that the HFHE technique can significantly improve the gas drainage efficiency by permeability enhancement. But the HFHE-enhanced gas drainage is unique to the other gas developments with high or low flow rates; a slight pressure gradient and short drainage time make the non-Darcy effect fail to generate cumulative influence on the HFHE-enhanced gas drainage.


Introduction
Due to the distinctive geological environment and the increasing mining level, most of China's coal mines in the main coal production bases face serious coal and gas outburst dangers. Degassing coal seams is an important method to mitigate this hazard and results in the beneficial recovery of a clean burning, low-carbon fuel resource. However, the permeability of China's coal seams is usually in the magnitude of 10 −4~1 0 −3 mD (except for Jincheng coalfield), which is four orders of magnitude lower than the US and three orders of magnitude lower than Australia [1,2]. The fairly low permeability challenges the feasibility, cost, and performance of coal mine methane (CMM) drainage.
Recently, hydraulic flushing hole enlargement technique has played an increasingly important role in solving the problem of poor efficiency associated with gas drainage, especially in soft coal seams [1,[3][4][5]. As shown in Figure 1(a), the HFHE-enhanced CMM drainage technique involves the use of high-pressure water jet in a borehole to flush the surrounding coal body. When an underground opening is excavated, the rock surrounding it must bear the load that was previously borne by the removed rock mass. If the rock is not strong enough, the borehole may be easily destroyed by the redistributed in situ stress (see Figure 1(b)), resulting in the in situ stress release and coal damage [6]. Various permeability experiments showed that during progressive damage to coal, permeability can increase by two to four orders of magnitude when a failure occurs [7,8]. The enhanced permeability of the surrounding coal is generally recognized as the main reason for the performance improvement of gas drainage.
Coal seams are typical dual porosity unconventional gas reservoirs that consist of coal matrix surrounded by intersecting fractures. Gas drainage would induce a gas migration in these dual porosity media: the gas diffusion in coal matrix is governed by Fick's law, and the viscous flow in coal fractures is normally described by Darcy's law [9]. Darcy's law describes a linear relationship between flow velocity and pressure gradient. However, when the flow velocity is greater than a certain level, the viscous flow will deviate from the linear relationship. This nonlinear relationship is defined as the non-Darcy flow [10,11].
The existence of the non-Darcy effect in the flow of fluids through porous media has been known for quite some time. The significance of non-Darcy flow has been emphasized in the literature associated with the good productivity of natural gas, coalbed methane, tight gas, and shale gas. Miskimins et al. [10] found that non-Darcy flow effects have an impact on the performance of a hydraulically fractured well even at low flow rates. Friedel and Voigt [12] investigated the effects of non-Darcy flow on fractured well productivity in tight-gas reservoirs and found that neglecting non-Darcy flow in welltest analysis can lead to an erroneous interpretation. Ye et al. [11] proposed a coupling model taking into account the coal permeability change and the variable β factor to study the non-Darcy flow behavior in coal seams. Barree and Conway [13] presented a new equation to describe the relationship between the rate and potential gradient for porous media flow over the entire range of the Reynolds number. Ramey [14] found that the non-Darcy effect could change the early pressure distribution around a gas well and suggested to consider this effect in estimating the variations in flow capacity in short-time gas well testing.
Though the non-Darcy effect in sandstone, limestone, and tight gas reservoir was extensively studied, few experimental and field studies were conducted on the non-Darcy flow in coal seams. Particularly, the HFHE technique can significantly enhance the coal permeability; however, non-Darcy flow effects in the HFHE-enhanced CMM drainage have been neglected in a majority of cases so far, without quantitative evaluations.
The objective of this paper is, therefore, to investigate the effects of non-Darcy flow on the HFHE-enhanced CMM productivity in a coal seam, by taking into account the gascoal interaction effect and Klinkenberg effect as reference objects. A fully coupled model of coal deformation and damage, non-Darcy flow of free methane gas in the fractures, and Fickian diffusion of adsorbed methane in the coal matrix was developed. A series of simulation scenarios were carried out to explore some key effects that impact the HFHEenhanced CMM productivity. Based on these efforts, the HFHE-induced improvement of gas drainage efficiency was analyzed, and the influences of the non-Darcy effect, gas-coal interaction effect, and Klinkenberg effect on the HFHE-enhanced gas drainage were evaluated and compared.

Non-Darcy Flow and Permeability-Dependent Non-Darcy Coefficient
The coal seam is a typical unconventional gas reservoir, which consists of a coal matrix surrounded by intersecting fractures. Gas migration is complicated in such a dualporosity system. This study simplifies the gas migration based on the conceptual model shown in Figure 2(a). The dual-porosity system has two different flow systems: the first phase involves the Fickian diffusion of adsorbed methane gas from the coal matrix into the free phase within voids in the coal seam, and the second phase is the viscous flow of free methane gas in the fracture network [9,15]. In this dual-porosity system, two pressures are present: one in the fractures, p f , and the other in the coal matrix blocks, p m . The flow in the fracture network is usually assumed to follow Darcy's law, which describes the velocity as a linear function of the pressure gradient: where k is the permeability, mD; μ is the dynamic methane viscosity, Pa·s; and v is the gas velocity, m/s. However, as shown in Figure 2(  2 Geofluids a certain range of velocities (and pressure gradients, correspondingly) for laminar flow, and for the velocities above or below some critical values, the flow becomes Forchheimer or pre-Laminar [16]. This nonlinear phenomenon is the so-called non-Darcy flow, which is primarily caused by the continuous deceleration and acceleration of fluid molecules traveling along a tortuous flow path through the interconnected pores and also in the proppant pack [17]. This non-Darcy flow can be described by adding a quadratic term of velocity to Darcy's law [11,18], which is often named as the Forchheimer equation as where ρ g is the gas density, kg/m 3 , and β is the non-Darcy flow coefficient, which is a reservoir-specific parameter and is correlating with permeability, m −1 . Cooke [19] proposed that β factor could be described as a power law of permeability based on the investigation of non-Darcy flow for fractures packed with multiple layers of sand: where a and b are the reservoir-specific parameters. Equation (3) makes it easier to estimate the non-Darcy effect by obtaining the reservoir permeability. Several studies have investigated the quantitative correlation between the β factor and the permeability of different reservoirs. As shown in Figure 3(a), the correlation for the experimental data from different sources and researchers is [12] β = 4 1 ⋅ 10 11 k −1 5 4 Permeability variation can significantly influence the non-Darcy effect. Unlike the conventional reservoir whose permeability variation is negligible, the coal permeability changes enormously during the enhanced CMM capture engineering. Primarily, the achievements of the enhanced CMM capture technologies are based on the evolution of the coal fracture network. As shown in Figure 3(b), experimental results show that during progressive deformation of coal, permeability first decreases as preexisting cleats close and then increases as new vertical dilatant microcracks are generated, which will suddenly increase by 3-4 orders of magnitude when a failure occurs [7]. At the same time, it is well known that coal permeability also changes due to the variation of effective stress and sorption-induced coal deformation during the production period. Thus, how the non-Darcy flow impacts the permeability-enhanced gas drainage should be evaluated, and a variable β factor based on the permeability model should be implemented.
Meanwhile, to incorporate the non-Darcy effect into numerical simulation more expediently, the Forchheimer equation can be reformulated as

Geofluids
This is the apparent Darcy's law, and a correction factor is implemented below to define the non-Darcy effect: where f T is the non-Darcy factor, which is a measure of the degree of non-Darcy flow through reservoirs [20]. If the β factor or velocity approaches zero, then the non-Darcy factor approaches one and the fluid flow obeys Darcy's law. k T is the equivalent permeability measured under non-Darcy flow condition (while k is the permeability measured under Darcy flow condition), mD, and it can be calculated by

A Mathematical Model for Enhanced CMM Capture
3.1. Gas Storage. In situ coal contains methane both in the coal matrix blocks as an adsorbed phase and as a free phase compressed in the fracture networks. Unlike conventional gas reservoirs, methane in coal is primarily stored as an adsorbed phase at near liquid densities [21], whose mass can be calculated by the popular Langmuir equation as where m a is the adsorbed gas mass per volume of coal, kg/m 3 ; V L is the Langmuir volume constant, kg/m 3 ; P L is the Langmuir pressure constant, MPa; ρ c is the density of coal, Figure 3: Non-Darcy flow coefficients for different reservoir rocks and coal permeability evolution during progressive deformation [7,12]. 4 Geofluids kg/m 3 ; and ρ s is the methane gas density at the standard condition, kg/m 3 and can be calculated by where M g is the molar mass of methane, kg/mol, and V M is the molar volume of methane at the standard condition, m 3 /mol. Free gas mass can be calculated by where m f is the free gas mass per volume of coal, kg/m 3 ; ϕ is the porosity of coal; and ρ f is the compressed gas density, kg/m 3 , and can be calculated by the ideal gas law: where R is the universal gas constant, J/(mol·K), and T is the temperature, K. The total gas mass m per volume of coal can be expressed as

Governing Equations for Gas Migration in a Coal
Seam. We will describe the gas migration in the coal seam by two balance equations: one for free gas in fractures and another for adsorbed gas in coal matrix blocks. By applying the two mass balance equations of gas in both the fracture network and the coal matrix, we have where t is the time, d, and Q m is the gas exchange rate per volume of coal, kg/(m 3 ·s). Gas desorbs in the coal block and then drains to the fractures governed by the methane concentration gradient; the gas exchange rate can be expressed by [22,23] where D is the gas diffusion coefficient, m 2 /s; σ c is the shape factor of coal matrix, m −2 ; and c m and c f are the concentrations of gas in the coal matrix and fractures, respectively, kg/m 3 , which can be derived by the ideal gas law: Moreover, for CMM reservoir modeling, σ c is related to sorption time τ, d. To facilitate computations, the sorption time is commonly used to approximate the diffusivity of coal matrix blocks. It has a reciprocal relationship with the diffusion coefficient and shape factor, and in general, it is numerically equivalent to the time during which 63.2% of the coal gas content is desorbed [24] Thus, the gas exchange rate Q m can be reformulated as Substituting (8), (9), and (18) into (14), we obtained the governing equation for the gas pressure change in the coal matrix blocks under the non-Darcy flow condition: Substituting (5), (10), (11), and (18) into (13), we obtained the governing equation for the gas pressure change in the fractures under the non-Darcy flow condition: Geomechanical Characterization of the HFHE-Enhanced CMM Drainage. Since avoiding rock mass failure is usually a major rock mechanics goal, knowledge of what happens after failure could seem to be of little interest. However, unlike the general rock engineering, the success of enhanced CMM capture engineering is based on the artificial failure behavior of coal. Based on a wide experience in the field of underground coal mining engineering and laboratory observations, the geomaterials are elastic-plastic materials and generally show strain-softening behavior after peak-load [25,26]. For the strain-softening behavior, the strength parameter gradually reduces from the peak state parameter to the residual state parameter after coal failure (see Figure 4). The transition from the peak to the residual state in the strain-softening zone is governed by the softening parameter. Based on internal variables, the most utilized softening parameter is a plastic shear strain, γ p , which is equal to 5 Geofluids the difference between the major and minor principal plastic strains: where ε p 1 and ε p 3 are the major and minor principal plastic strains, respectively.
Thus, the strength parameters can be defined by the bilinear functions of plastic shear strain [27]: where ω represents one of the strength parameters; ω e and ω r represent one the strength parameters at the elastic zone and residual zone, respectively; and γ p r is the critical plastic shear strain from which the residual behavior starts.
The Mohr-Coulomb (M-C) criterion is one of the most widely used failure criteria in geomechanical engineering due to its mathematical simplicity and clear physical meaning. However, the M-C criterion plot in the pi-plane is an irregular hexagon with three "sharp" corners, as opposed to a smooth function, e.g., Drucker-Prager (D-P) failure criterion, which limits the numerical implementation of M-C criterion [28]. To eliminate this problem, the usual strategy applied in numerical simulation is to make circumscribed and inscribed combinations of the Coulomb and Drucker-Prager criteria [29]. The linkage between the M-C and the D-P criteria is achieved by computing the external or circumscribed D-P model with the Coulomb parameters.
The M-C matching D-failure criterion is given by [28,30] where J 2 is the second invariant of the stress deviator tensor; I 1 is the first invariant of the stress deviator tensor; and λ and κ are the material constants, which can be computed as follows: where φ is the internal friction angle of the coal,°, and c is the cohesion of the coal, MPa.

Evolution of Coal Porosity and Permeability.
Permeability is the key parameter in determining gas migration in the coal. It is well known that coal permeability is governed by a range of fracture characteristics which are sensitive to the stress state of the coal [31]. To enhance the coal permeability, we need to make the stress redistributed in the coal. Thus, during the progressive deformation of coal, the fracture characteristics can change due to normal stress-induced closures or openings and shear stress-induced dilations, then inducing the coal permeability evolution [32]. For simplicity, the permeability varying with stress state during the progressive deformation can be defined by trilinear functions of stress and plastic shear strain (see Figure 4). At the elastic zone, the permeability evolution with the stress change can be expressed as an exponential decline [33][34][35]: where k 0 represents the permeability at the initial stress state, mD; c f is the fracture volume compressibility, MPa −1 ; and Θ represents the volumetric stress, MPa. As deviatoric stress increases, new fractures, favorably oriented along the direction of the maximum principal stress, will be created, which will change the permeability from decline to increase [7,36]. As shown in Figure 4, when a coal failure occurs after the peak stress is reached, the permeability suddenly increases by several orders of magnitude. The permeability increases linearly with the softening deformation after peak stress and keeps constant at the residual state [37,38]. Thus, the transition of permeability from the peak stress to the residual state in the strain-softening zone and the residual zone can also be defined by the softening parameter as follows: where ζ represents the jump coefficient of permeability (JCP), which can be calculated by k b /k a , where k b and k a represent the permeability values measured at the peak state and residual state, respectively (see Figure 4). Due to the decrease of gas pressures in both the fractures and the coal matrix with time during the CMM production,  6 Geofluids the permeability will vary with the competing effects of effective stress and sorption-induced coal deformation. The evolution of permeability during CMM production can be described by the modified Palmer-Mansoori (PM) model [39]: where k * ∞0 represents the initial absolute permeability before CMM production, mD; ς represents the Klinkenberg factor, Pa, where ς = 0 251k −0 36 ∞ ; f p f , p m , ε s represents the competing effects of effective stress and sorption-induced coal deformation on coal permeability; ε s represents the sorption-induced strain of coal; ε L is the Langmuir volumetric strain at infinite pore pressure; K is the bulk modulus of coal, MPa, where Young's modulus of coal, MPa; υ is Poisson's ratio of coal; and χ is the effective stress coefficient. Based on the widely used cubic law, the porosity change during the CMM production can be defined by [40] By combining (25), (26), and (27) In this section, the model developed above was applied to quantify the non-Darcy effect on the performance of a drainage borehole. A series of simulation scenarios were conducted and solved with COMSOL Multiphysics, a commercial software for the implementation of partial differential equations (PDEs). As listed in Table 1, six study cases were carried out to study four factors influencing gas drainage performance, including the stress redistribution effect, gas-coal interaction effect, Klinkenberg effect, and non-Darcy effect. The HFHE-enhanced permeability and β factor characterization were achieved, and their significances on gas drainage were discussed. Viscous flow velocity and non-Darcy factor distributions were studied, and their relationship with the non-Darcy effect was explored. Finally, the non-Darcy effect in the HFHEenhanced CMM drainage was evaluated by comparing to the gas-coal interaction effect and Klinkenberg effect. Particularly, the present work only focused on single-phase flow, and the water influence was not involved at this stage. The governing equations of diffusion and non-Darcy flow were solved by the PDE module, and the deformation equations for the HFHE-enhanced permeability distribution were solved by the solid mechanics module. The geometry and boundary conditions of the simulation model are illustrated in Figure 5. A square region with a length of 20 m which was simplified from the HFHEenhanced CMM drainage engineering was set as the geometry model. The length of the geometry model was implemented to ensure it was longer than the typical borehole spacing in Yangquan coalfield. The drainage borehole was placed at the corner, whose radius was set based on the study cases according to Table 1. Suitable boundary conditions must be applied to the simulation model. For the solid deformation model, symmetry boundary conditions were applied to the left and bottom sides, and the stress conditions were applied to the top and right sides. The initial displacement in the domain was zero. In the gas migration model, a constant gas pressure of 87 kPa was applied to the borehole, and no flow boundary conditions were applied to   Table 2.

HFHE-Enhanced Permeability Characterization.
To better understand the impact of the HFHE-induced permeability evolution on the non-Darcy effect in coal seams, a series of sensitivity studies on the key parameter (i.e., jump coefficient of permeability) have been conducted. When a borehole is excavated, the stress redistribution near wellbore occurs causing stress changes around the wellbore compared to the in situ or far-field stress. Once the redistributed stress exceeds the coal failure strength, the elastic-plastic secondary stress distribution will be formed and the surrounding coal of a borehole will generate damage. Figure 6 illustrates damage variable distribution around a borehole under different conditions. The radius of a widely used borehole in an underground coal mine in China is 0.054 m, which can be enlarged to 0.5 m by the HFHE technique. According to (26), it can be concluded that the generation of coal damage and redistribution of in situ stress will significantly influence the permeability distribution around a borehole. The damage zone of the borehole whose radius is 0.5 m is greater than that whose radius is 0.054 m, indicating that the HFHE technique can remarkably improve the permeability of the surrounding coal. Figure 7 illustrates the relative permeability distribution around a borehole under different conditions. The relative permeability is equal to the result of the redistributed permeability divided by its initial value (i.e., k/k 0 ). The relative permeability data in a square with 3 m long was analyzed. For a convenient comparison, the color ranges for different conditions keep the same, which ranges from 0.84 to 2600. It is obvious that the closer to the borehole wall, the greater 9 Geofluids the permeability, and the permeability will keep constant beyond the borehole effect region. Under the same geomechanical condition, it is the borehole radius that controls the scale of the borehole effect region. The effect region of the borehole whose radius is 0.5 m is quite bigger than the one whose radius is 0.054 m, indicating that the HFHE technique can significantly enlarge the affect region of in situ stress and permeability.
The permeability and non-Darcy flow coefficient (β factor) along the detection line under all the 6 conditions were extracted to further analysis on the influences of borehole radius and JCP on permeability (see Figure 8). Figures 8(a) and 8(b) represent the result of boreholes whose radiuses are 0.054 m and 0.5 m, respectively. It can be found that the permeability and β factor show obvious partitioning feature, and the inflection point is closely related to the exterior boundary of the borehole effect region. When the borehole radius increases from 0.054 m to 0.5 m, the radius of the borehole effect region increases from 0.2 m to 1.8 m. The permeability gradually decreases with the increasing distance to the borehole wall, and the decay rate gradually increases. The permeability in most part of the borehole effect region is greater than the initial value, which is almost equal to the initial value beyond the borehole effect region. Under the same geomechanical condition, the permeability in the borehole effect region is mainly controlled by JCP. The higher the JCP, the greater the permeability. JCP is closely related to the structural characteristics of coal and can be determined by the seepage characteristics during the complete stress-strain compression test of coal. In this study, not a specified JCP was provided; three JCP values covering the potential permeability change were analyzed.
β factor is a reservoir-specific parameter which is a power law function of permeability with a negative index (see (4)), and thus, the β factor change has an opposite trend to the permeability change (see Figure 8). The β factor gradually increases with the increasing distance to the borehole wall and then tends to go back to the initial value. Due to the negative power function relation between β factor and permeability, the β factor can maximally increase by four orders of magnitude when the permeability only increases hundreds of times (JCP = 100). β factor is one of the three main impacts that control the equivalent permeability k T , and the HFHE technique can significantly influence the β factor. Thus, whether the non-Darcy flow effect in the HFHE-enhanced CMM drainage can be neglected or not should be seriously evaluated.

HFHE-Enhanced Gas Flow
Characterizations. The discussion above in Section 4.2 focused on the HFHEenhanced permeability characterizations. In this section, the HFHE-enhanced gas flow characterizations will be discussed based on the comparison between case A and case B. Particularly, as the borehole effect region of case A is far smaller compared to the solution domain, the redistributed permeability in this small region will lead to singularity, and thus, the HFHE-induced permeability redistribution was involved in case A. Figure 9 presents the pressure distributions of case A and case B at different drainage times. The color range in Figure 9 represents the gas pressure ranges from 87 kPa (average drainage subpressure) to 2.5 MPa (initial gas pressure). The minimum gas pressures are those surrounding the drainage borehole, and the gas pressure gradually increases with the increasing distance to the borehole wall. The pressure decay region induced by borehole expands with the increase of drainage time. It is obvious that the pressure decay region of case B is bigger than that of case A when borehole drainage works. After drainage of 100 d, the pressure decay region of case A was still small. In contrast, due to the HFHEenhanced permeability in the borehole effect region, a relatively bigger pressure decay region was formed within drainage of 1 d. The quick formation of the pressure decay region will lead to plenty of gas emission from the borehole in a short time. In the field, a large amount of gas will be  10 Geofluids released from the borehole in a very short time due to the HFHE technique. Thus, a gas collector is needed for the HFHE-enhanced gas drainage in the field. The numerical simulation results are consistent with the gas emission behaviors of the HFHE-enhanced CMM drainage. Figure 10 shows the distributions of gas pressure along the detection line and relative gas content in the detection zone, respectively. In consideration of the distance between two adjacent boreholes in the HFHE-enhanced gas drainage, the relative gas content in a square with 5 m long is analyzed. The gas content in the detection zone can be calculated by the discrete surface integral method.
As shown in Figure 10(a), it is obvious that the gas pressure distribution of case B is different from that of case A. For case B, there is a distinct pressure decay part with 1.5 m long, which is below 0.4 MPa after drainage of 1 d. The gas pressure can decrease to the average drainage subpressure in this part after drainage of 100 d. In the borehole drainage affect region, the gas pressure for case B is much lower than that of case A.
The difference in drainage efficiency between case A and case B can be more easily observed from the relative gas content evolution with time (see Figure 10(b)). The gas stored in the coal seam is gradually released with gas drainage. The relative gas content for case B first quickly decreases and then gradually decreases. At any time, the decay rate of gas content for case B is higher than that of case A. Therefore, the gas stored in the coal seam decreases only about 2.5% for case A and 22.5% for case B after drainage of 100 d, indicating that the HFHE technique can significantly improve the drainage efficiency.
The HFHE-induced improvement of gas drainage efficiency is benefited from two aspects. One can be directly associated with the HFHE-induced permeability 11 Geofluids enhancement. The other reason is that the quick decrease of gas pressure in the permeability enhancement zone will form an "equivalent borehole" with a larger radius. It has already been proved that the increment in the borehole radius can benefit the borehole drainage efficiency.

Non-Darcy Effect, Gas-Coal Interaction Effect, and
Klinkenberg Effect. Figure 11 presents the gas pressure distributions along the detection line of different case studies at various times. Figure 12 presents the relative gas content evolution with time in the detection zone of different case studies. The detection line and detection zone are set the same as those in Section 4.3.1. Based on the comparison of gas pressure and gas content of different case studies, the influences of non-Darcy effect, gas-coal interaction effect, and Klinkenberg effect on the HFHE-enhanced gas drainage can be evaluated.
The gas-coal interaction effect can be analyzed by comparing the results of case B and case C. It can be found that during the whole drainage period (100 d), the gas pressure is higher when the gas-coal interaction effect is taken into account. The coal permeability decreases with the gas pressure decrease due to the gas-coal interaction effect; i.e., the effective stress is dominant compared with the sorptioninduced deformation in the drainage period (100 d). After drainage of 100 d, the difference of relative gas content between case B and case C is not more than 0.5%.
The non-Darcy effect was implemented in case E. At the very initial stage (t = 0.05 d), the gas pressure distribution of case E is almost the same to that of case C, indicating that the impact of the non-Darcy effect is of the same order of magnitude as that of the gas-coal interaction effect. But it is difficult to distinguish the difference in the gas pressure distribution between case B and case E except for the very initial stage (t = 0.05 d). The relative gas content evolution with time of case E is almost equal to that of case B for the whole drainage period. Thus, it can be concluded that the non-Darcy effect only has a slight influence on the HFHEenhanced gas drainage. However, hundreds of papers have been published revealing that the non-Darcy effect can significantly influence the good productivity of natural gas, coalbed methane, tight gas, and shale gas [10][11][12]41]. Compared to these gas developments, the non-Darcy effect only has a slight influence on the HFHE-enhanced gas drainage performance, and the main reasons are as follows.
Though the coal permeability can increase by hundreds or thousands of times through the HFHE technique, the permeability-enhanced area is still limited. At the same time, the initial gas pressure is only 2.5 MPa, which is far lower than that of tight-gas reservoirs, shale reservoirs, and etc. Thus, the gas pressure in the permeability-enhanced area decreases very quickly. The slight pressure gradient cannot maintain a high gas velocity. That is why the non-Darcy effect has been only observed at the very initial stage. Moreover, the initial coal permeability is too small to generate non-Darcy flow in the coal seam beyond the enhanced area.
Miskimins et al. [10] found that the non-Darcy effect could reduce cumulative gas production by up to 18.1% over a ten-year period in low permeability reservoirs. However, due to the engineering characteristic, the HFHE-enhanced gas drainage usually sustains a short period (less than 1 year). Thus, the non-Darcy effect can not generate cumulative influence on the HFHE-enhanced gas drainage.
In contrast, the Klinkenberg effect has a significant influence on the HFHE-enhanced gas drainage. The gas pressure of case D is lower than that of case B in the influenced area of a borehole for the whole drainage period. The Klinkenberg effect may have a significant impact on gas flow behavior, especially in low permeability media [42]. Though the gas pressure differences between case B and case D become smaller with the increase of drainage, time increases in the permeability-enhanced region. The Klinkenberg effect does impact the other parts of the borehole influenced region. Thus, after drainage of 100 d, the difference of relative gas content between case B and case D can reach about 4%. Case F involves all the three effects (i.e., non-Darcy effect, gas-coal interaction effect, and Klinkenberg effect), and at this condition, the gas pressure distribution and gas content evolution are similar to those of case D, indicating that the Klinkenberg effect is dominant.

Conclusions
This study numerically investigated the effects of non-Darcy flow on the HFHE-enhanced CMM productivity in the coal seam through a fully coupled numerical model. The gas migration in the coal seam involves two flow systems: the non-Darcy flow of free methane gas in the fracture network and the Fickian diffusion of adsorbed methane gas from the coal matrix into the free phase within fractures. Based on the analysis of geomechanical characterization of the HFHE technique, a permeability model for describing the stress redistribution and gas-coal interaction-induced permeability evolution has been proposed. By using this numerical model, a series of study cases were carried out to explore some key effects that impact the HFHEenhanced CMM productivity.
Based on the results obtained, the primary findings reached in this study are as follows: (1) The HFHE technique can significantly enlarge the affected region of in situ stress and permeability. When the borehole radius increases from 0.054 m to 0.5 m, the radius of the borehole effect region increases from 0.2 m to 1.8 m. The permeability in most part of the borehole effect region is greater than the initial value, which is mainly controlled by JCP.
The β factor can maximally increase by four orders of magnitude when the permeability only increases hundreds of times (2) The quick formation of the pressure decay region induced by the HFHE technique will lead to plenty of gas emission from the borehole in a short time, which is in accordance with field performance. The gas stored in the coal seam is gradually released due to gas drainage, and the HFHE technique can increase the gas content decay ratio from 2.5% to 22.5% after drainage of 100 d, indicating that the HFHE technique can significantly improve the CMM drainage efficiency (3) The HFHE-induced improvement of gas drainage efficiency is benefited from two aspects. One can be directly associated with the HFHE-induced permeability enhancement. The other reason is that the quick decrease of gas pressure in the permeability enhancement zone will form an "equivalent borehole" with a larger radius (4) The influences of the non-Darcy effect, gas-coal interaction effect, and Klinkenberg effect on the HFHE-enhanced gas drainage have been evaluated and compared. Without considering the gas-coal interaction effect and the Klinkenberg effect, the calculated drainage performance will be underestimated by 0.5% and overestimated by 4% after drainage of 100 d, respectively. However, a slight pressure gradient and short drainage time make the non-Darcy effect fail to generate cumulative influence on the HFHE-enhanced gas drainage. Comparing to the gas-coal interaction effect and Klinkenberg effect, the non-Darcy effect only has a slight influence on the HFHE-enhanced gas drainage and can be ignored Simulation of the non-Darcy effect in a given producing situation is not simple, and the numerical work presented in this study is not intended to be a substitute for 13 Geofluids detailed analysis in a specific situation. However, it does provide a tool that can be used for evaluating the potential non-Darcy effect in a given HFHE-enhanced CMM production situation.

Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.

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