Experimental Study on the Permeability Evolution and Gas Flow Law of Post-Strength Soft Coal

: Permeability is an essential indicator for predicting gas drainage yield and preventing mine gas disasters, which is signiﬁcantly inﬂuenced by the stress paths and the integrity of coal. Conventional research on permeability mainly focused on the permeability evolution of initial undamaged or fractured (prefabricated fractures) coal under various stress paths; little attention has been paid to post-strength coal (stress-induced damage), especially for soft coal. To determine the permeability evolution and gas ﬂow law of post-strength soft coal samples under various stress paths, we used the experimental method combined with the numerical method in this study. The results showed that when the conﬁning pressure and axial pressure of post-strength soft coal samples were unloaded, the permeability increased by 1.25–1.32 times; when the coal samples were loaded into the secondary damage, the permeability ﬁrst decreased and then increased. The simulation part in this study found that the development of the fracture of coal samples under triaxial compression was divided into four stages. Gas ﬂow law of post-strength soft coal was signiﬁcantly inﬂuenced by fracture locations, and the gas pressure and gas ﬂow ﬁeld near the fracture were disturbed.


Introduction
China is one of the countries with the most severe gas disasters in the world [1]. According to the statistics, from 2013 to 2017, there were 1945 coal mine accidents in China, of which 218 were caused by gas, accounting for 11.21% of the total [2]. Coal is a porous medium with abundant internal pores [3]. To recognize the gas flow law in coal, many theories have been developed, e.g., Darcy's permeation, Fick's diffusion, the diffusion-permeation mixture, etc. [4][5][6][7]. Based on the above theories, researchers found that the flow law of gas in coal is mainly determined by the gas pressure and the permeability of the coal seam; the permeability is significantly influenced by the stress paths and the integrity of coal [8,9]. Thus, it is necessary to investigate the permeability evolution of different integrity (undamaged, fractured, post-strength) states of coal under various stress paths (loading, cyclic loading, unloading).
Previous studies on the permeability changes with different stress paths and different integrity states of coal are mainly divided into four aspects. First, the permeability change rules of initial

Apparatus
A triaxial creep-seepage test system for coal and rock equipment (Henan Polytechnic University, Jiaozuo, China) was employed to conduct the triaxial seepage tests. The instrument consisted of a press and model system, a gas supply system, a pressure and temperature test system, an outlet metering system, a computer acquisition and control processing system, etc. This equipment had high stiffness (5000 kN/mm), high load-sensing accuracy (2.5 N), and high deformation-sensing accuracy (0.75 mm). Its axial load could be controlled by force or displacement. Relevant parameters, including stress, strain, deformation, and gas pressure, were collected and recorded automatically during the triaxial seepage tests. The test device is shown in Figure 2. (a)

Apparatus
A triaxial creep-seepage test system for coal and rock equipment (Henan Polytechnic University, Jiaozuo, China) was employed to conduct the triaxial seepage tests. The instrument consisted of a press and model system, a gas supply system, a pressure and temperature test system, an outlet metering system, a computer acquisition and control processing system, etc. This equipment had high stiffness (5000 kN/mm), high load-sensing accuracy (2.5 N), and high deformation-sensing accuracy (0.75 mm). Its axial load could be controlled by force or displacement. Relevant parameters, including stress, strain, deformation, and gas pressure, were collected and recorded automatically during the triaxial seepage tests. The test device is shown in Figure 2.

Apparatus
A triaxial creep-seepage test system for coal and rock equipment (Henan Polytechnic University, Jiaozuo, China) was employed to conduct the triaxial seepage tests. The instrument consisted of a press and model system, a gas supply system, a pressure and temperature test system, an outlet metering system, a computer acquisition and control processing system, etc. This equipment had high stiffness (5000 kN/mm), high load-sensing accuracy (2.5 N), and high deformation-sensing accuracy (0.75 mm). Its axial load could be controlled by force or displacement. Relevant parameters, including stress, strain, deformation, and gas pressure, were collected and recorded automatically during the triaxial seepage tests. The test device is shown in Figure 2. (a)

Experimental Principle and Scheme
Suppose there is a rock with a cross-sectional area A and a length L, which is clamped in the core holder, as shown in Figure 3. The flow rate Q is measured by passing a fluid having a viscosity of μ through the center of the rock at a pressure difference ΔP. The experiment proves that the volume flow rate Q through the core per unit time is proportional to the pressure difference ΔP core crosssectional area A and inversely proportional to the core length L and the viscosity of the fluid μ: (1) Equation (1) is the so-called "Darcy equation." It can be seen from the Equation (1) that A and L are the geometrical dimensions of the rock and ΔP is the external condition. When the external conditions, square sizes, and fluid properties are constant, the fluid throughput Q depends on the proportionality constant K. We refer to K as the permeability of the rock, then Equation (1) can be rewritten as: The previous discussion was based on the incompressible fluid (liquid). The gas permeability we designed was based on the compressed fluid (liquid), so the "Darcy equation" needs to be modified to apply it.
It is well known that the characteristics of compressible gas are that when the pressure increases, the fluid can be compressed; when the pressure is lowered, the fluid expands; when the temperature is constant, the expansion of the fluid obeys Boyle's law. If considering the simplest flat line seepage, the inlet pressure is set as P1, and the outlet pressure is set as P2. When the pressure changes from P1 to P2, the volume of the gas necessarily changes, so the flow rate also varies. Therefore, it is necessary to consider substituting the average volumetric flow rate into the Darcy equation.
If the gas expansion is regarded as an isothermal process, according to Boyle's law: Q= 2P 0 Q 0 P 1 +P 2 (5) Figure 2. Three-axis creep-seepage test system for coal and rock. (a) Image of the triaxial cell; (b) schematic diagram of the three-axis creep-seepage test system.

Experimental Principle and Scheme
Suppose there is a rock with a cross-sectional area A and a length L, which is clamped in the core holder, as shown in Figure 3. The flow rate Q is measured by passing a fluid having a viscosity of µ through the center of the rock at a pressure difference ∆P. The experiment proves that the volume flow rate Q through the core per unit time is proportional to the pressure difference ∆P core cross-sectional area A and inversely proportional to the core length L and the viscosity of the fluid µ: Equation (1) is the so-called "Darcy equation." It can be seen from the Equation (1) that A and L are the geometrical dimensions of the rock and ∆P is the external condition. When the external conditions, square sizes, and fluid properties are constant, the fluid throughput Q depends on the proportionality constant K. We refer to K as the permeability of the rock, then Equation (1) can be rewritten as: Appl. Sci. 2020, 10, x FOR PEER REVIEW 6 of 28 where: P is the mean pressure, MPa; Q is the average volume flow under average pressure, mL/s; P0 is the atmospheric pressure, MPa; Q0 is the gas flow at atmospheric pressure, mL/s. From the above analysis, the correction of the Darcy formula of the compressible fluid is: where μg is the viscosity of the gas, Pa s. To simulate the permeability evolution with repeated increase and decrease of the stress on the post-strength coal seam in the process of coal mining, the test scheme adopted in this test was as follows: When the gas pressure was 0.7 MPa (this value was based on the geological conditions of the Da-ping Coal Mine), the confining pressure was 2 MPa, 3 MPa, and 4 MPa, and the axial stress was loaded until the coal samples were damaged. Then, the permeability evolution in the process of unloading axial stress and confining pressure and the permeability evolution in the process of loading axial stress to the secondary damage were measured. To make the experiment method clearer, we list the experimental scheme in Table 1. It should be noted that to determine the The previous discussion was based on the incompressible fluid (liquid). The gas permeability we designed was based on the compressed fluid (liquid), so the "Darcy equation" needs to be modified to apply it.
It is well known that the characteristics of compressible gas are that when the pressure increases, the fluid can be compressed; when the pressure is lowered, the fluid expands; when the temperature is constant, the expansion of the fluid obeys Boyle's law. If considering the simplest flat line seepage, the inlet pressure is set as P 1 , and the outlet pressure is set as P 2 . When the pressure changes from P 1 to P 2 , the volume of the gas necessarily changes, so the flow rate also varies. Therefore, it is necessary to consider substituting the average volumetric flow rate into the Darcy equation.
If the gas expansion is regarded as an isothermal process, according to Boyle's law: where: P is the mean pressure, MPa; Q is the average volume flow under average pressure, mL/s; P 0 is the atmospheric pressure, MPa; Q 0 is the gas flow at atmospheric pressure, mL/s. From the above analysis, the correction of the Darcy formula of the compressible fluid is: where µ g is the viscosity of the gas, Pa s. To simulate the permeability evolution with repeated increase and decrease of the stress on the post-strength coal seam in the process of coal mining, the test scheme adopted in this test was as follows: When the gas pressure was 0.7 MPa (this value was based on the geological conditions of the Da-ping Coal Mine), the confining pressure was 2 MPa, 3 MPa, and 4 MPa, and the axial stress was loaded until the coal samples were damaged. Then, the permeability evolution in the process of unloading axial stress and confining pressure and the permeability evolution in the process of loading axial stress to the secondary damage were measured. To make the experiment method clearer, we list the experimental scheme in Table 1. It should be noted that to determine the permeability of coal samples, we had to measure the gas flow rate. The error of permeability was mainly affected by the flowmeter, and the flowmeter we adopted in this test was 0.5 L with an accuracy of ±10 mL/min, then we could obtain the absolute deviation of permeability as ±0.045 mD. The data were collected automatically every 6 s during the experiment; due to the large amount of data, the data interval we used for analyzing was 18 s. Note: C i is the initial confining pressure; A i is the peak stress; A u is the axial stress after unloading; C u is the confining pressure after unloading; C sl is the confining pressure after the secondary loading; A sl is the axial stress after the secondary loading.

Coal Sample M1
Coal sample M1 was under the experimental condition in which the gas pressure was 0.7 MPa and the confining pressure 2 MPa. The axial stress was loaded until the coal sample M1 was damaged. Then, the permeability evolution in the process of unloading the axial stress and confining pressure and  It can be seen from Figure 4a that after the coal sample M1 was damaged, the axial stress was unloaded from 12 MPa to 1 MPa at a rate of 20 N/S. As the axial stress decreased, the permeability gradually increased. When the axial stress was unloaded to 1 MPa, the permeability increased to about 1.08 mD, which was about 1.3 times that before unloading the axial stress.
It can be seen from Figure 4b that when the confining pressure was unloaded from 2 MPa to 1 MPa, the permeability gradually increased. When the confining pressure was unloaded to 1 MPa, the permeability increased to about 1.14 mD, which was approximately 1.32 times that before unloading.
It can be seen from Figure 4c,d that as the axial stress increased, the permeability gradually decreased. When the axial stress increased to a certain extent, the coal sample M1 underwent secondary damage, and the permeability began to increase steadily. The permeability changed to approximately U-shaped during the secondary damage process. The permeability of the coal sample M1 after initial damage was about 0.85 mD, and the permeability minimum value in the secondary damage process was about 1.05 mD. It can be seen that the permeability of the coal sample increased as a whole during the secondary damage process.

Coal Sample M2
Coal sample M2 was under the experimental condition in which the gas pressure was 0.7 MPa and the confining pressure 3 MPa. The axial stress was loaded until the coal sample M2 was damaged. Then, the permeability evolution in the process of unloading the axial stress and confining pressure and loading the axial stress to the secondary damage of the coal sample M2 was measured. The following curves were obtained, as shown in Figure 5. It can be seen from Figure 4a that after the coal sample M1 was damaged, the axial stress was unloaded from 12 MPa to 1 MPa at a rate of 20 N/S. As the axial stress decreased, the permeability gradually increased. When the axial stress was unloaded to 1 MPa, the permeability increased to about 1.08 mD, which was about 1.3 times that before unloading the axial stress.
It can be seen from Figure 4b that when the confining pressure was unloaded from 2 MPa to 1 MPa, the permeability gradually increased. When the confining pressure was unloaded to 1 MPa, the permeability increased to about 1.14 mD, which was approximately 1.32 times that before unloading.
It can be seen from Figure 4c,d that as the axial stress increased, the permeability gradually decreased. When the axial stress increased to a certain extent, the coal sample M1 underwent secondary damage, and the permeability began to increase steadily. The permeability changed to approximately U-shaped during the secondary damage process. The permeability of the coal sample M1 after initial damage was about 0.85 mD, and the permeability minimum value in the secondary damage process was about 1.05 mD. It can be seen that the permeability of the coal sample increased as a whole during the secondary damage process.

Coal Sample M2
Coal sample M2 was under the experimental condition in which the gas pressure was 0.7 MPa and the confining pressure 3 MPa. The axial stress was loaded until the coal sample M2 was damaged. Then, the permeability evolution in the process of unloading the axial stress and confining pressure and loading the axial stress to the secondary damage of the coal sample M2 was measured. The following curves were obtained, as shown in Figure 5.  It can be seen from Figure 5a that after the coal sample M2 was damaged, the permeability increased when the axial stress was unloaded from 15 MPa to 1 MPa, and the increased rate of permeability was slow first and then fast. The permeability was lower than that of coal sample M1, which was under the confining pressure of 2 MPa. When the axial stress was reduced to 1 MPa, the permeability increased to about 0.9 mD, which was about 1.25 times that before unloading the axial stress.
It can be seen from Figure 5b that when the confining pressure was unloaded from 3 MPa to 1 MPa, the radial restraint force of the coal sample decreased, and the permeability range increased approximately linearly. When the confining pressure was unloaded to 1 MPa, the permeability increased to about 1.15 mD, which was about 1.32 times that before unloading the confining pressure.
It can be seen from Figure 5c that when the confining pressure of the M2 coal sample was loaded from 1 MPa to 2 MPa, the fractures of the coal sample M2 were closed again, and its permeability was approximately linearly reduced. When the confining pressure was increased by 1 MPa, the permeability was reduced by about 1.17 times.
It can be seen from Figure 5d that when the axial stress was loaded from 1 MPa to 9 MPa, the coal sample underwent secondary damage. The peak intensity of the primary damage of the coal sample M2 was 2.33 times that of the secondary damage. During the secondary damage process, the change trend of permeability decreased first and then increased. Compared with Figure 5a, it can be seen that the permeability of the coal sample M2 after initial damage was about 0.72 mD, and the permeability minimum value in the secondary damage process was about 0.68 mD. The average permeability in the secondary damage process was higher than that during the initial damage process.

Coal Sample M3
Coal sample M3 was under the experimental condition in which the gas pressure was 0.7 MPa and the confining pressure 4 MPa. The axial stress was loaded until the coal sample M3 was damaged. It can be seen from Figure 5a that after the coal sample M2 was damaged, the permeability increased when the axial stress was unloaded from 15 MPa to 1 MPa, and the increased rate of permeability was slow first and then fast. The permeability was lower than that of coal sample M1, which was under the confining pressure of 2 MPa. When the axial stress was reduced to 1 MPa, the permeability increased to about 0.9 mD, which was about 1.25 times that before unloading the axial stress.
It can be seen from Figure 5b that when the confining pressure was unloaded from 3 MPa to 1 MPa, the radial restraint force of the coal sample decreased, and the permeability range increased approximately linearly. When the confining pressure was unloaded to 1 MPa, the permeability increased to about 1.15 mD, which was about 1.32 times that before unloading the confining pressure.
It can be seen from Figure 5c that when the confining pressure of the M2 coal sample was loaded from 1 MPa to 2 MPa, the fractures of the coal sample M2 were closed again, and its permeability was approximately linearly reduced. When the confining pressure was increased by 1 MPa, the permeability was reduced by about 1.17 times.
It can be seen from Figure 5d that when the axial stress was loaded from 1 MPa to 9 MPa, the coal sample underwent secondary damage. The peak intensity of the primary damage of the coal sample M2 was 2.33 times that of the secondary damage. During the secondary damage process, the change trend of permeability decreased first and then increased. Compared with Figure 5a, it can be seen that the permeability of the coal sample M2 after initial damage was about 0.72 mD, and the permeability minimum value in the secondary damage process was about 0.68 mD. The average permeability in the secondary damage process was higher than that during the initial damage process.

Coal Sample M3
Coal sample M3 was under the experimental condition in which the gas pressure was 0.7 MPa and the confining pressure 4 MPa. The axial stress was loaded until the coal sample M3 was damaged.
Then, the permeability evolution in the process of unloading the axial stress and confining pressure and loading the axial stress to the secondary damage of the coal sample M3 was measured. The following curves were obtained, as shown in Figure 6.
Appl. Sci. 2020, 10, x FOR PEER REVIEW 9 of 28 Then, the permeability evolution in the process of unloading the axial stress and confining pressure and loading the axial stress to the secondary damage of the coal sample M3 was measured. The following curves were obtained, as shown in Figure 6. It can be seen from Figure 6a that the axial stress of the coal sample M3 was loaded to 20 MPa under the condition of 4 MPa confining pressure, and then, the initial damage occurred. In the process of unloading axial stress from 20 MPa to 1 MPa, the permeability of coal sample M3 increased slowly. When the axial stress was unloaded to 1 MPa, the permeability increased by about 1.3 times that before unloading the axial stress.
According to Figure 6b, when the confining pressure was unloaded from 4 MPa to 1 MPa, the fractures binding force of coal sample M3 decreased, and the permeability increased gradually. When the confining pressure was unloaded to 1 MPa, the permeability was about 1.25 times that before unloading the confining pressure.
It can be seen from Figure 6c that when the confining pressure of the coal sample M3 was loaded from 1 MPa to 3 MPa, the fractures of the coal sample M3 were re-closed, and the gas seepage channel became smaller. The permeability range decreased approximately linearly, and the permeability decreased about 1.22 times when the confining pressure was increased by 2 MPa.
It can be seen from Figure 6d that when the axial stress was loaded from 1 MPa to 13 MPa, the coal sample M3 underwent secondary damage. The peak intensity of the initial damage of the coal sample M3 was 1.76 times that of the secondary damage. During the whole process, the permeability of the coal sample M3 first decreased and then increased. Compared with Figure 6a, it can be seen that the permeability of the coal sample M2 after initial damage was about 0.64 mD, and the permeability minimum value in the secondary damage process was about 0.62 mD. However, in general, the permeability of coal sample M3 in the secondary damage process was higher than that during the initial damage process, which was consistent with the permeability evolution of coal samples M1 and M2 in the same processes. It can be seen from Figure 6a that the axial stress of the coal sample M3 was loaded to 20 MPa under the condition of 4 MPa confining pressure, and then, the initial damage occurred. In the process of unloading axial stress from 20 MPa to 1 MPa, the permeability of coal sample M3 increased slowly. When the axial stress was unloaded to 1 MPa, the permeability increased by about 1.3 times that before unloading the axial stress.
According to Figure 6b, when the confining pressure was unloaded from 4 MPa to 1 MPa, the fractures binding force of coal sample M3 decreased, and the permeability increased gradually. When the confining pressure was unloaded to 1 MPa, the permeability was about 1.25 times that before unloading the confining pressure.
It can be seen from Figure 6c that when the confining pressure of the coal sample M3 was loaded from 1 MPa to 3 MPa, the fractures of the coal sample M3 were re-closed, and the gas seepage channel became smaller. The permeability range decreased approximately linearly, and the permeability decreased about 1.22 times when the confining pressure was increased by 2 MPa.
It can be seen from Figure 6d that when the axial stress was loaded from 1 MPa to 13 MPa, the coal sample M3 underwent secondary damage. The peak intensity of the initial damage of the coal sample M3 was 1.76 times that of the secondary damage. During the whole process, the permeability of the coal sample M3 first decreased and then increased. Compared with Figure 6a, it can be seen that the permeability of the coal sample M2 after initial damage was about 0.64 mD, and the permeability minimum value in the secondary damage process was about 0.62 mD. However, in general, the permeability of coal sample M3 in the secondary damage process was higher than that during the initial damage process, which was consistent with the permeability evolution of coal samples M1 and M2 in the same processes.

Permeability Evolution
The permeability evolution of coal samples with time is shown in Figure 7.

Permeability Evolution
The permeability evolution of coal samples with time is shown in Figure 7.  The table could be obtained b()y the statistics of the test data of each coal sample, as shown in Table 2. As shown in Figure 7 and Table 2, the peak stress of coal samples M1, M2, and M3 was 14 MPa, 21 MPa, and 23 MPa, and the corresponding permeability of the peak stress point was 0.84 mD, 0.72 mD, and 0.64 mD. The permeability evolution trend of M2 and M3 was almost the same, but was different from M1; this was because M1 had no secondary loading confining pressure stage. In the secondary damage period, the rate of permeability evolution was significantly higher than in other conditions. The maximum value and the minimum value of permeability decreased with the growth of confining pressure, such as at Point E, the permeability of M1, M2, and M3 was 1.15 mD, 0.72 mD, and 0.63 mD, which showed a decreasing trend.
As shown in Figure 8, when the axial stress was unloaded, the permeability of coal samples increased, and the permeability decreased with the growth of confining pressure. When the axial stress was unloaded to 1 Mpa, the permeability of M1, M2, and M3 was 1.06 mD, 0.89 mD, and 0.80 mD. When the confining pressure was unloaded, the permeability increased; when confining pressure was unloaded to 1 MPa, the permeability of M1, M2, and M3 was 1.14 mD, 1.12 mD, and 1.03 mD. It can be seen that the permeability of M1 was higher than that of M2 and M3 under the same confining pressure. When the confining pressure was loaded, the permeability decreased; when the confining pressure was loaded to 2 MPa, the permeability of M2 and M3 was 0.95 mD and 1.1 mD. When the coal samples were loaded to the secondary damage, the permeability first decreased, then increased, and the permeability still decreased with the increase of confining pressure. From the above analysis, we could know under any condition that the permeability was M1 > M2 > M3, and this was due to that fact that as the confining pressure increased, the micropores in the coal body were compacted, resulting in a decrease in porosity. The permeability evolution with axial stress and confining pressure in the unloading process showed an exponential relation; however, in the loading confining pressure process, the relationship was linear.

Model Establishment
The coal samples' damage processes were also the cumulative damage and fracture process of the internal fracture generation, expansion, and penetration of the coal matrix. Studying the

Model Establishment
The coal samples' damage processes were also the cumulative damage and fracture process of the internal fracture generation, expansion, and penetration of the coal matrix. Studying the mechanism of fractures propagation and penetration is essential for studying the evolution of permeability. The numerical simulation model was established in the RFPA (Realistic Failure Process Analysis) 2D with a size of W50 mm × H100 mm, which was divided into 120 × 240 grids. The Mohr-Coulomb strength criterion was selected. The loading mode was stress loading with each step loaded with 0.1 MPa, and the bottom was fixed both vertically and horizontally. The fracture development mechanism of the coal samples under different confining pressure was compared and studied. The established model is shown in Figure 9, and the basic parameters of the model are shown in Table 3.

Simulation Results and Analysis
The axial stress was loaded at the initial stress of 10 MPa and 0.1 MPa per step, until the coal samples were damaged under the confining pressures of 2 MPa, 3 Mpa, and 4 MPa, respectively. The fractures development law in the damaged process of coal samples was obtained. The maximum shear stress distribution and acoustic emission of coal samples under different confining pressure are shown in Figure 10.

Simulation Results and Analysis
The axial stress was loaded at the initial stress of 10 MPa and 0.1 MPa per step, until the coal samples were damaged under the confining pressures of 2 MPa, 3 Mpa, and 4 MPa, respectively. The fractures development law in the damaged process of coal samples was obtained. The maximum shear stress distribution and acoustic emission of coal samples under different confining pressure are shown in Figure 10. As shown in Figure 10a1, when the confining pressure of the coal sample was 2 MPa, the axial stress was loaded from 10 MPa, 0.1 MPa was loaded per step, and a total of 37 steps was loaded. When the model was loaded to Step 10, the coal matrix produced distinct damaged units. When the model was loaded to Step 25, the damaged units were distributed in the coal matrix. When the model was loaded to Step 37, the stress concentration at the fracture tip led to the expansion of the fracture; the damaged units passed through each other and eventually aggregated to form a macroscopic fracture; thus, the coal sample was unstable and damaged.
As shown in Figure 10a2, the diameter of the circle represents the relative size of the acoustic emission intensity. The white is the acoustic emission generated by compression shear damage, and the red is the acoustic emission generated by tensile damage. At the initial stage of loading, acoustic emission was generated randomly in the coal matrix and mainly was in compression shear failure. With the increase of load, acoustic emission caused by compression shear failure was distributed throughout the whole coal sample. When the model was loaded to Steps 37-36, the acoustic emission of tensile failure began to appear, then a large number of acoustic emission counts of tensile failure occurred until the damage of the coal sample occurred; the acoustic emission gathered at the fractures and finally formed the acoustic emission band consistent with the fracture surface. As shown in Figure 10a1, when the confining pressure of the coal sample was 2 MPa, the axial stress was loaded from 10 MPa, 0.1 MPa was loaded per step, and a total of 37 steps was loaded. When the model was loaded to Step 10, the coal matrix produced distinct damaged units. When the model was loaded to Step 25, the damaged units were distributed in the coal matrix. When the model was loaded to Step 37, the stress concentration at the fracture tip led to the expansion of the fracture; the damaged units passed through each other and eventually aggregated to form a macroscopic fracture; thus, the coal sample was unstable and damaged.
As shown in Figure 10a2, the diameter of the circle represents the relative size of the acoustic emission intensity. The white is the acoustic emission generated by compression shear damage, and the red is the acoustic emission generated by tensile damage. At the initial stage of loading, acoustic emission was generated randomly in the coal matrix and mainly was in compression shear failure. With the increase of load, acoustic emission caused by compression shear failure was distributed throughout the whole coal sample. When the model was loaded to Steps 37-36, the acoustic emission of tensile failure began to appear, then a large number of acoustic emission counts of tensile failure occurred until the damage of the coal sample occurred; the acoustic emission gathered at the fractures and finally formed the acoustic emission band consistent with the fracture surface.
The fractures produced in this process would have a significant impact on gas permeability. When the model was loaded to the Steps 37-36, the damaged units in the coal sample were relatively distinct, and the coal sample was in the plastic stage. At this time, due to the increase of the porosity of the coal sample, the gas permeability would show an increasing trend. With the continuous growth of the load, the porosity of the coal sample would continue to increase, and the corresponding gas permeability would gradually increase until the model was loaded to Steps 37-44, and the fractures formed a large fluid channel; thus the gas permeability reached the peak value.
As shown in Figure 10b1, when the confining pressure of the coal sample was 3 MPa, the axial stress was loaded from 10 MPa, 0.1 MPa was loaded per step, and a total of 76 steps were loaded. The damage process was similar to that of σ = 2 MPa. It is described as follows: When the model was loaded to Steps 25, the coal matrix produced significantly damaged units. When the model was loaded to Steps 60, the damaged units were distributed in the coal matrix. When the model was loaded to Steps 76, the stress concentration at the fracture's tip led to the expansion of the fracture, and the damaged units passed through each other and eventually aggregated to form macroscopic fractures; thus, the coal sample was unstable and damaged.
As shown in Figure 10b2, when the model was loaded to Steps 76-27, the lower-left corner of the model produced evident acoustic emission, and the lower-left corner presented micro-fractures. Due to the stress concentration at the end of the fractures, the fractures expanded upward along the lower-left corner as the load continued to increase. When the model was loaded to Steps 76-31, the acoustic emission penetrated the entire model at 45 degrees substantially along the lower-left corner of the fractures. When the model was loaded to Steps 76-35, the model was damaged, and the acoustic emission was concentrated at the fractures, thus forming a visible shear-damage acoustic emission band.
Similarly, in this process, the gas permeability also changed with the changes in the porosity of the coal sample. When the model was loaded to Steps 76-27, the damaged units in the coal sample converged to form tiny fractures at the lower-left corner. With the continuous increase of the load, the damaged units in the coal matrix continued to unite and connect, and the porosity of the coal sample continued to increase; thus, the gas permeability gradually increased. When the model was loaded to Steps 76-35, the gas permeability reached the peak due to the formation of the through fractures.
As shown in Figure 10c1, when the confining pressure of the coal sample was 4 MPa, the axial stress was still loaded from 10 MPa, 0.1 MPa was loaded per step, and a total of 104 steps were loaded. The damaged process is described as follows: When the model was loaded to Steps 35, the coal matrix produced distinct damaged units. When the model was loaded to Steps 60, the damaged units were distributed in the coal matrix. When the model was loaded to Steps 104, significant fractures occurred at the lower end of the coal sample. The coal sample underwent damage at Steps 104-20, producing an "X"-type fracture.
As shown in Figure 10c2, the fractures first occurred in the lower part of the coal sample. When the model was loaded to Steps 62, this part began to produce a red acoustic emission circle. When the model was loaded to Steps 104, the acoustic emission was concentrated in this part, and the model was unstable and damaged; the acoustic emission finally formed an acoustic emission band along the fracture surface.
In this process, the gas permeability changed with the changes in the porosity of the coal sample. When the model was loaded to Steps 104-12, the damaged units in the coal sample converged in the middle and lower part of the coal matrix, and the gas permeability began to increase. When the model was loaded to Steps 104-16, the damaged area was enlarged, and the gas permeability continued to grow. When the model was loaded to Steps 104-20, the coal sample was unstable and damaged. At this time, the corresponding gas permeability reached a peak value.
The AE (acoustic emission) energy of coal samples during loading is shown in Figure 11, and the AE (acoustic emission) counts are shown in Table 4. 1.2x10 -9 Energy(MJ) Step Energy (a) 1.2x10 -9 1.4x10 -9 1.6x10 -9 Energy(MJ) Step Energy (b) Step Energy (c)   It can be seen from Figure 11 and Table 4 that because of the confining pressure, the primary fractures inside the coal sample were closed by different degrees of compaction, which improved the integrity and rigidity of the coal sample. The existence of confining pressure led to the small acoustic emission energy in the initial loading. Compared with the energy change of 2 MPa, 3 MPa, and 4 MPa, it can be seen that in the initial loading stage, the acoustic emission energy decreased with the increases of the confining pressure. For example, when the confining pressure was 4 MPa, there was almost no acoustic emission energy. With the rise of the axial stress, the acoustic emission energy also showed an upward trend, because there was a small scale of new fracture growth and generation in the coal sample. When the coal sample was in the elastic-plastic stage, there were more new fractures in the coal sample; at the same time, a large amount of strain energy was accumulated in the coal sample, and the acoustic emission activity was further strengthened. When the coal sample entered the stage of fracture development, in addition to the expansion of the original fractures, a large number of new fractures was produced. These fractures converged and linked up, resulting in the formation of the macro fracture surface, and the interaction between fractures was intensified, for example Step 37 of σ = 2 MPa, Step 76 of σ = 3 MPa, and Step 104 of σ = 4 MPa were all stages of the macro-destructive energy's sudden increase. When the coal sample entered the residual deformation stage, the acoustic emission signal returned to a lower level. At this stage, the damage of the coal sample was completed, and the fractures were closed under the action of pressure.
Combining the development of the fracture with acoustic emission during the damage process of coal samples, we divided the development of the fracture into four stages: (1) Micro-damage nucleation and growth stage characterized by the micro-damages being randomly distributed in the coal sample; the scales were not much different from each other, and the statistics were independent of each other.

Study of the Gas Flow Law in Fractured Coal Samples
Although only the data of three coal samples were cited above, a plurality of experiments was carried out due to the influence of the experimental success rate. Most of the damaged coal samples obtained through experiments had a shear failure through the fractures, and the shear fracture locations would be different. Some were located in the upper part of the coal sample, some in the middle, and also some in the lower part of the coal sample. Therefore, it was necessary to analyze the distribution of the gas pressure and flow field at different fracture positions of the coal samples by numerical simulation and to obtain the inherent law through comparative analysis.

Established Model
Gas flow models with varying locations of the fracture of coal samples were established in the COMSOL numerical simulation software. The model had a width of 50 mm and a height of 100 mm. Gas flowed in from the upper boundary of the model and flowed out from the lower boundary; the upper boundary was set as the pressure inlet, the pressure being 0.7 MPa; the lower boundary was set as the pressure outlet, the pressure being 0.1 MPa. The boundary on both sides was set as no flow; the fracture boundary selected fracture flow and defined parameters such as the porosity and permeability of fracture. The established models and relevant parameters are shown in Figure 12 and Table 5.
upper boundary was set as the pressure inlet, the pressure being 0.7 MPa; the lower boundary was set as the pressure outlet, the pressure being 0.1 MPa. The boundary on both sides was set as no flow; the fracture boundary selected fracture flow and defined parameters such as the porosity and permeability of fracture. The established models and relevant parameters are shown in Figure 12 and Table 5.

Upper Fracture Model
It can be seen from Figure 13 that the internal pressure of the model changed due to the generation of the upper fracture. On both sides of the fracture, the pressure on the left side was generally higher than that on the right side, especially at the fracture exit. The fracture made the pressure isoline of the model distribute in the fracture exit as wave-shaped bending. The closer to the fracture exit, the higher the curvature of the isoline bending. The reason for the above phenomenon was the upper fracture on the left side, which induced the gas flow on the left side of the model and made the gas pressure on the left side of the model higher.

Upper Fracture Model
It can be seen from Figure 13 that the internal pressure of the model changed due to the generation of the upper fracture. On both sides of the fracture, the pressure on the left side was generally higher than that on the right side, especially at the fracture exit. The fracture made the pressure isoline of the model distribute in the fracture exit as wave-shaped bending. The closer to the fracture exit, the higher the curvature of the isoline bending. The reason for the above phenomenon was the upper fracture on the left side, which induced the gas flow on the left side of the model and made the gas pressure on the left side of the model higher.

Upper Fracture Model
It can be seen from Figure 13 that the internal pressure of the model changed due to the generation of the upper fracture. On both sides of the fracture, the pressure on the left side was generally higher than that on the right side, especially at the fracture exit. The fracture made the pressure isoline of the model distribute in the fracture exit as wave-shaped bending. The closer to the fracture exit, the higher the curvature of the isoline bending. The reason for the above phenomenon was the upper fracture on the left side, which induced the gas flow on the left side of the model and made the gas pressure on the left side of the model higher.  It can be seen from Figure 14 that as the gas flowed mostly along the fracture, the flow had a concentrated effect, resulting in a significant decrease in the gas flow on the left and right sides of the fracture, and almost no flow in the upper left portion of the fracture. When the gas flowed down the fracture to the exit position, the flow was gradually distributed evenly to fill the entire model. At the interface between the fracture and the matrix, the gas flowed vertically into the matrix from the fracture surface.
Appl. Sci. 2020, 10, x FOR PEER REVIEW 21 of 28 It can be seen from Figure 14 that as the gas flowed mostly along the fracture, the flow had a concentrated effect, resulting in a significant decrease in the gas flow on the left and right sides of the fracture, and almost no flow in the upper left portion of the fracture. When the gas flowed down the fracture to the exit position, the flow was gradually distributed evenly to fill the entire model. At the interface between the fracture and the matrix, the gas flowed vertically into the matrix from the fracture surface.

Middle Fracture Model
According to the analysis of Figure 15, the pressure distribution decreased from top to bottom. On both sides of the middle fracture, the pressure on the left side was higher than that on the right side, especially at the fracture outlet. The fracture made the pressure isoline inlet and outlet of the model curve wave shaped. The closer to the inlet and outlet of the fracture, the higher the curvature of the contour curve was, and the curvature of the outlet pressure contour was higher than that of the inlet pressure contour. The influence range of the middle fracture on the pressure contour was about 0.06 m along the Y direction, and the influence rate of the outlet was higher than that of the inlet.

Middle Fracture Model
According to the analysis of Figure 15, the pressure distribution decreased from top to bottom. On both sides of the middle fracture, the pressure on the left side was higher than that on the right side, especially at the fracture outlet. The fracture made the pressure isoline inlet and outlet of the model curve wave shaped. The closer to the inlet and outlet of the fracture, the higher the curvature of the contour curve was, and the curvature of the outlet pressure contour was higher than that of the inlet pressure contour. The influence range of the middle fracture on the pressure contour was about 0.06 m along the Y direction, and the influence rate of the outlet was higher than that of the inlet. As shown in Figure 16, the gas flow velocity at the inlet and outlet of the fracture was much faster than that at other locations. The upper-right velocity distribution and the lower-left velocity distribution were symmetrical, and the gas velocity distribution of the two parts was more significant. The upper left velocity distribution and the lower right velocity distribution were balanced, and the gas velocity distribution of the two parts was smaller. It can be seen that, due to the effect of fractures, most of the gas flowed down along the upper left of the model through the fracture channel, and a small part of the gas flowed vertically from the upper right to the lower right. After the gas flowed to the lower part of the fracture for a certain distance, with the weakening of the fracture effect, it gradually started to flow evenly. As shown in Figure 16, the gas flow velocity at the inlet and outlet of the fracture was much faster than that at other locations. The upper-right velocity distribution and the lower-left velocity distribution were symmetrical, and the gas velocity distribution of the two parts was more significant. The upper left velocity distribution and the lower right velocity distribution were balanced, and the gas velocity distribution of the two parts was smaller. It can be seen that, due to the effect of fractures, most of the gas flowed down along the upper left of the model through the fracture channel, and a small part of the gas flowed vertically from the upper right to the lower right. After the gas flowed to the lower part of the fracture for a certain distance, with the weakening of the fracture effect, it gradually started to flow evenly. As shown in Figure 16, the gas flow velocity at the inlet and outlet of the fracture was much faster than that at other locations. The upper-right velocity distribution and the lower-left velocity distribution were symmetrical, and the gas velocity distribution of the two parts was more significant. The upper left velocity distribution and the lower right velocity distribution were balanced, and the gas velocity distribution of the two parts was smaller. It can be seen that, due to the effect of fractures, most of the gas flowed down along the upper left of the model through the fracture channel, and a small part of the gas flowed vertically from the upper right to the lower right. After the gas flowed to the lower part of the fracture for a certain distance, with the weakening of the fracture effect, it gradually started to flow evenly. As shown in Figure 17, the gas pressure decreased from top to bottom. In the lower part of the model, nearly half of the space in the model was in a low-pressure state due to the penetration of the lower fracture and the bottom boundary. The contour of the pressure distribution at the upper end of the fracture was wave shaped, and the curvature was smaller as it went up. The pressure on the

Lower Fracture Model
As shown in Figure 17, the gas pressure decreased from top to bottom. In the lower part of the model, nearly half of the space in the model was in a low-pressure state due to the penetration of the lower fracture and the bottom boundary. The contour of the pressure distribution at the upper end of the fracture was wave shaped, and the curvature was smaller as it went up. The pressure on the left side of the fracture was higher than that on the right side, and the influence of the lower fracture on the pressure contour was about 0.04 m in the Y direction.
(c) (d) Figure 16. Distribution of gas velocity in the middle fracture model. (a) Gas velocity in the coal matrix; (b) gas velocity field in the coal matrix; (c) gas velocity in middle fracture; (d) gas velocity filed in middle fracture.

Lower Fracture Model
As shown in Figure 17, the gas pressure decreased from top to bottom. In the lower part of the model, nearly half of the space in the model was in a low-pressure state due to the penetration of the lower fracture and the bottom boundary. The contour of the pressure distribution at the upper end of the fracture was wave shaped, and the curvature was smaller as it went up. The pressure on the left side of the fracture was higher than that on the right side, and the influence of the lower fracture on the pressure contour was about 0.04 m in the Y direction. As shown in Figure 18, the gas flow in the upper half of the model was relatively uniform as a whole, and the flow velocity at the fracture inlet was relatively high. The flow velocity in the upper half of the model was higher than that in the lower half. Part of the gas entered the fracture channel from the fracture inlet, and part of the gas entered the fracture channel from the matrix perpendicular to the fracture surface. The velocity field was uniform in the upper half, pointing vertically downward. The velocity field began to skew towards the entrance at a position of about 0.02 m above the fracture entrance. Due to the existence of the lower fracture, most of the gas flowed out of the model along the fracture channel, and only a small part flowed out from the right lower side of the model. As shown in Figure 18, the gas flow in the upper half of the model was relatively uniform as a whole, and the flow velocity at the fracture inlet was relatively high. The flow velocity in the upper half of the model was higher than that in the lower half. Part of the gas entered the fracture channel from the fracture inlet, and part of the gas entered the fracture channel from the matrix perpendicular to the fracture surface. The velocity field was uniform in the upper half, pointing vertically downward. The velocity field began to skew towards the entrance at a position of about 0.02 m above the fracture entrance. Due to the existence of the lower fracture, most of the gas flowed out of the model along the fracture channel, and only a small part flowed out from the right lower side of the model.     Figure 20a shows that the velocity of the fracture outlet of the upper fracture model was about 0.0039 m/s, which was much higher than that at other locations along the transversal line. The pressure of the fracture outlet of the upper fracture model was about 690 KPa, which was due to the connection between the fracture and the upper boundary of the model, resulting in the downward movement of the high-pressure area. The pressure at the fracture outlet was higher than that at other  Figure 20a shows that the velocity of the fracture outlet of the upper fracture model was about 0.0039 m/s, which was much higher than that at other locations along the transversal line. The pressure of the fracture outlet of the upper fracture model was about 690 KPa, which was due to the connection between the fracture and the upper boundary of the model, resulting in the downward movement of the high-pressure area. The pressure at the fracture outlet was higher than that at other positions at the same level. Figure 20b,c show that the velocity of the fracture inlet and outlet of the middle fracture model was about 0.0028 m/s and 0.0029 m/s, respectively. The velocity along the transversal line showed that the closer to the fracture inlet and outlet, the higher the flow velocity. The pressure of the fracture inlet of the middle fracture model was about 470 KPa, which was lower than that at other locations along the transversal line; the pressure of the fracture outlet of the middle fracture model was about 450 KPa, which was higher than other locations along the transversal line. Figure 20d shows that the velocity of the fracture inlet of the lower fracture model was about 0.0052 m/s, which was higher than that at other locations along the transversal line; the pressure was about 50 KPa, which was lower than that at other locations along the transversal line.

Conclusions
We used the experimental methods combined with the numerical method to investigate the permeability evolution of post-strength soft coal. First, the coal samples with different confining pressures (2 MPa, 3 MPa, 4 MPa) were loaded until they reached peak strength and were damaged. Then, the permeability evolution of post-strength soft coal samples during unloading axial stress and confining pressure and during loading to secondary damage was studied. Lastly, we adopted numerical simulation software to acquire the development mechanism of the internal fracture network of coal samples and the gas flow law of coal samples with various fracture locations. The main conclusions can be summarized as below: (1) When the confining pressure and axial pressure of post-strength coal samples were unloaded, the permeability increased by 1.25-1.32 times, and the permeability evolution with axial stress and confining pressure in the unloading process showed an exponential relation; however, in the loading confining pressure process, the relationship was linear. (2) When the coal samples were loaded into the secondary damage process, the change trend of permeability first decreased and then increased, and the permeability during the secondary damage of the coal samples was higher than during the initial damage.