Numerical Simulation of Coupled Thermal-Hydrological- Mechanical-Chemical Processes in the Spontaneous Combustion of Underground Coal Seams

In this study, we develop a fully coupled thermal-hydrological-mechanical-chemical (THMC) model to analyze the spontaneous combustion process of underground coal seams, focusing on investigating the influences of the pressure difference between oxygen and coal, the rate of coal-oxygen reaction heat, and the activation energy. The simulation results show that as oxygen propagates into the coal seams, the coal-oxygen reaction causes the spontaneous combustion of coal to heat. The consumption of oxygen leads to an increase in oxygen consumption along the way and a decrease in gas pressure. The permeability near the right boundary increases while significantly reducing the area far away from the right boundary as the predominant effect of spontaneous combustion. Additionally, a sensitivity study shows that a more considerable pressure difference and coal-oxygen reaction heat contribute to promoting the coal temperature, while the activation energy has a slight effect. Moreover, an increase in coal-oxygen reaction heat and activation energy accelerates the oxygen consumption rate and thus causes a lower oxygen concentration. Overall, the results provide a basis for the prediction and prevention of coal seam spontaneous combustion.


Introduction
The problem of coal seam spontaneous combustion is recognized as an urgent problem in coal production, which has the potential possibility leading to production accidents such as gas explosions, as well as severe economic and environmental problems [1][2][3][4][5][6]. In China, nearly 100-200 million tons of coal are lost to spontaneous combustion every year. Additionally, CO 2 content produced by spontaneous combustion accounts for 2-3% of global carbon emissions [7].
The process of coal spontaneous combustion is believed to be a complex physical and chemical process, including the contact and flow of coal and oxygen in the cracks of the coal seam, the accumulation and diffusion of heat when oxygen is sufficient, the distribution of temperature fields, and the transfer of energy [5,6,8]. Wu et al. [9] studied the influence of air and dust on the com-bustion process near the spontaneous combustion point of coal by numerical analysis of the kinetic equation. Su et al. [10] analyzed the influence of coal seam inclination and airflow on coal spontaneous combustion by establishing a three-dimensional model. Xia et al. [11] showed a coupled modeling of coal seam heat flow solidification. The effect of gas expansion changes in pressure gradient on gas seepage and heat transfer is investigated. Zongqing [12] studied the competitive adsorption relationship among multicomponent gases in the spontaneous combustion process of underground coal seams and the influence of gas concentration on the development of coal pores. Xia et al. [13] studied the interaction of multiphysics coupling during the spontaneous combustion of methane-rich coal seams. The simulation results showed that the resolved methane gas dilutes the oxygen concentration and hinders oxygen migration in the porous medium. Liu and Qin [14] studied the mechanism of spontaneous combustion in the coal seam from experiments, which demonstrates that the spontaneous combustion phenomenon is the interaction and coupling of the internal pressure, velocity, oxygen concentration field, gas, and coal solid temperature field in the coal seam. Hongfen et al. [15] considered coupling factors of the mechanical, temperature, seepage, mass transfer, and other multiphysical fields. They summarized the location and change law of both high temperature and high displacement areas in the process of coal spontaneous combustion.
In the combustion process, the rock mechanics simulation of underground spontaneous coal seams shows that the shrinkage of coal seams and changes in mechanical properties lead to the collapse of the overlying rock layer and larger fissures accelerating the combustion process [3,16]. Based on the gas percolation and temperature diffusion in spontaneous combustion coal seams, Huang et al. [17] concluded that the coal seam permeability is greater than 10 -9 m 2 in order for the gas to fully diffuse and circulate into the coal spontaneous combustion process. Meanwhile, the heat transfer method is dominated by heat convection in shallow coal seams. Heat transfer is dominated by heat conduction in the combustion of deep coal seams. Wen [18] established a super large-scale spontaneous combustion experiment platform, which considers the actual on-site situation and the critical value of the coal temperature change rate during the combustion process. Song and Kuenzer [19] considered the limited influence of oxygen convection and diffusion and mass transfer process on the reaction rate in coal fire oxidation and combustion and deduced the calculation formula of coal oxidation and combustion rate in the hightemperature stage of coal fire. Fusheng et al. [20] studied the influence of the pore structure characteristics of coal samples on the spontaneous combustion of coal by an experimental method. Wessling et al. [21] separated the spontaneous combustion process of coal into different stages. They performed a multiphysics coupling numerical simulation and investigated the influence of coal seam permeability anisotropy on coal spontaneous combustion. Wen et al. [22] used different gas flow equations to construct three coupling models. The influence of gas seepage on the temperature distribution in the natural process of coal is explored. Dudzińska and Cygankiewicz [23] analyzed the impact of the gas produced and the gas adsorbed by the coal body through experiments. Lu [24] established a mathematical model to investigate the influence of plastic deformation on porosity and combustion effect caused by coal spontaneous combustion, which is used to analyze the on-site examples.
Considering the coupling of multiphysics is an attractive and effective method for studying coal seam spontaneous combustion. However, due to the complexity of the theoretical equations, the differences in boundary conditions and initial conditions significantly cause variations in the process of each physics. In this paper, we establish a fully coupled thermal-hydrological-mechanical-chemical (THMC) model to study the spontaneous combustion process of coal seams. Then, we investigate the influences of the pressure difference between oxygen and coal, the rate of coal-oxygen reaction heat, and the activation energy.

Governing Equations of Coupled THMC Model
The spontaneous combustion process of underground coal seams includes geomechanical deformation, gas seepage, oxygen component transmission, and heat energy transmission. The changes in temperature cause the deformation of coal and thus result in porosity and permeability evolution. The variation in porosity and permeability affects the behavior of flow and heat. This is a complicated two-way coupling relationship. The proposed model in this paper is derived based on the following assumptions: (1) Coal is a homogeneous and isotropic elastic medium (2) The deformation of coal satisfies the assumption of small deformation (3) The influence of gas and moisture in coal on spontaneous combustion of coal is neglected (4) The gas inside the pores and fractures of the coal is an ideal gas. The viscosity of gas keeps constant during the whole process

Coal Deformation Equations.
Based on the constitutive relations of poroelasticity, the governing equation for coal deformation that considers desorption-induced strain and temperature-induced strain can be expressed as follows [25,26]: where σ ′ and σ is effective and total stress, p is pore pressure, and C is the tangential stiffness matrix. The relationship between total strain ε and displacement is shown as The thermal strain and desorption-induced volumetric strain are calculated by where ΔT is the temperature change and α T is the thermal expansion coefficient; ε L and p L represent Langmuir strain constant and Langmuir pressure constant, respectively; b 1 and b 2 are the pressure and temperature coefficients, respectively; and T 0 is the reference temperature.

Gas Flow Equations.
In the dual-porosity and single permeability (DPSP) model, the mass conservation law of the 2 Geofluids coal matrix system can be described as the following: where m m represents the methane content per unit volume of the matrix, φ f is the porosity of fracture, and Q is the mass exchange between coal matrix and fracture, which can be defined as where χ is the shape factor and M g is molar mass of gas.
In this study, we embed the new proposed nonlinear diffusion model into Equation (10), which depends on diffusion time, gas pressure, and temperature. The total methane content in the coal matrix consists of two components: free gas and adsorbed gas. Thus, the value of m m can be given as follows: where φ m is the porosity of matrix, ρ is gas density, ρ a is the gas density under standard condition, ρ c is coal density, and V sg represent the adsorbed gas content in coal matrix, which is defined by Langmuir volume equation: where V L is the Langmuir volume constant.
Gas transport in fracture obeys Darcy's law and mass conservation law, so the governing equation for gas seepage in the fracture is expressed as [27] where k is cleat permeability, μ is the velocity of gas, and   3 Geofluids Q is the source or sink of gas. In the model, the value of Q equals to the diffusion amount of gas in the matrix. Additionally, coal deformation has a significant impact on fracture porosity and permeability. Based on Reference [28] and cubic law, the porosity and permeability varying with mean effective stress are shown as Equations (9) and (10), respectively, where σ′ is effective stress and φ f0 and k 0 are the initial porosity of fracture and permeability, respectively. Substituting Equations (9) and (10) into Equation (8), we obtain where S = ðφ f − α f Þ exp ð−Δσ ′ /KÞ.

Oxygen Transfer
Equations. The mass conservation of oxygen component flow in porous media can be expressed as [11,13] where D k is the gas diffusion coefficient, C O 2 are oxygen concentration and its reference value, A is preexponential factors, and E a is the activation energy.

Heat Transfer Equations.
Considering thermal dilatation of gas and coal, thermal convection, thermal conduction, and gas adsorption energy, the governing equation for heat transfer can be described as [29] ∂ ρC ð Þ M T Â Ã ∂t where K g is the bulk modulus of gas, α g ðα g = 1/TÞ is the thermal expansion coefficient of gas, Q h is oxidation reaction heat, ε v is volumetric strain, q st is the isosteric heat of adsorption, C g is gas-specific heat constant, λ s and λ g are thermal conductivities of coal and gas, respectively, and ðρCÞ M represents specific heat capacity of gas-filled coal, which can be expressed as where C s is the coal-specific heat constant.

Model Setup
In the proposed model, there are five variables ðp m , p f , c k , T, ε v , or σ ′ Þ for the coupled mathematic model. All the variables are interactive, and all the governing equations are nonlinear partial differential equations (PDE). COMSOL Multiphysics software [31] is adopted to solve the complex problem, which provides a powerful PDE-based modeling environment. In this study, the coal deformation is calculated by the solid mechanics module with Equation (1). Four general PDE modules are utilized to address gas dynamic diffusion, gas flow, concentration transfer, and heat transfer with Equations (11), (12), and (15), respectively. We use the time-stepping method of the implicit backward differentiation formula (BDF), which is an implicit solver that adopts backward differentiation formulas with variable discretization order and automatic step-size selection. In the time step, the Newton-like method is used to solve the algebraic equations with automatic linearization. Then, the resulting linear equations are solved employing multifrontal massively parallel sparse (MUMPS) solvers. The geometry of the two cases and mesh discretization of the simulation domain are shown in Figures 1 and 2, respectively. The size of the simulation domain is 5 m × 5 m. The simulation area is discretized into 6282 triangle elements and 200 edge elements. The initial temperature, pressure, and oxygen concentration are 300 K, 1 atm, and  Geofluids 0 mol/m 3 . The right boundary is assigned as the oxygen injection boundary with 9.375 mol/m 3 . The pressure and temperature of the right-side boundary are 2 MPa and 300 K, while the left-side boundary has the same pressure and temperature values as the initial conditions. The remaining boundaries are imposed with thermal insulation and no-flow conditions. The domain has boundary loadings of 6 MPa on the left and upper boundaries for the mechanical boundaries. The lower and right boundaries are constrained with zero displacements. The numerical simulation is carried out in two steps. In the first step, we calculated the initial equilibrium state with a ramped loading under the in situ stress, pressure, and temperature conditions. In the second step, the above simulated results are obtained as the initial conditions for the second step. The loads and boundary conditions on a transient model are ramped up from values that are consistent with the initial values. The smoothed step functions are adopted to ramp the boundary condition. The specific simulation parameters are shown in Table 1. Figure 3 shows the spatial distributions of the temperature at t = 1, 5, 10, and 20 days. The simulated result shows that temperature increases with the continuous oxygen propagation into the coal seams. Moreover, the range of oxygen influence expands with the simulation time, which is about 1 meter far away from the right boundary after 20 days. Figure 4 shows the distributions of gas pressure at t = 1 and 20 days. The coal-oxygen reaction causes the spontaneous combustion of coal to heat up and oxygen consumption along the way. The heat-up of coal spontaneous combustion further aggravates the coal-oxygen reaction, resulting in an increase in oxygen consumption along the way and a decrease in gas pressure. The hightemperature area is mainly concentrated in the area near the right boundary. Figure 5 shows the distributions of permeability ratio at t = 1, 5, 10, and 20 days. The expansion of the coal matrix caused by spontaneous combustion of coal causes the reduction of coal fracture opening. In contrast, the shrinkage of the coal matrix caused by gas dissociation results in the increase of the opening of coal fracture. The enhancement in permeability near the right boundary is observed, while the permeability in the area far away from the right boundary significantly reduces. Figure 6 shows the evolution of temperature, oxygen concentration, and permeability ratio at monitoring point A. After 20 days of spontaneous combustion of coal, the temperature generally increases up to about 360 K. The oxygen concentration at point A decreases with time after a short decrease. The maximum oxygen concentration inside the coal is about 3.65 kJ/mol 3 . The change of coal permeability is controlled by the combination of spontaneous combustion and shrinkage of the coal matrix. The permeability ratio reduced continually as the result of the predominated effect of spontaneous combustion of coal.

Effect of Pressure Difference.
In this section, we investigate the influence of pressure difference between oxygen and coal. We perform three cases in which the right boundary is imposed with 2 MPa, 3 MPa, and 4 MPa, respectively. Figure 7 shows the evaluation of temperature and oxygen concentration with different pressure differences. It can be clearly seen that the greater the pressure difference causes, the higher the air leakage rate and greater oxidation heating rate of coal. The maximum temperature at point A increases from around 360 K to 480 K.

Geofluids
Case c has the largest oxygen concentration of 0.72 kJ/mol as more oxygen propagates inside the coal for the higher flow or transfer velocity.

Effect of Oxidation Reaction
Heat. In this section, we carry out three cases of oxidation reaction heat Q n = 300 kJ/mol, 350 kJ/mol, and 400 kJ/mol to the analysis of effect of heat released by coal oxidation. Figure 8 shows the evaluation of temperature and oxygen concentration with different values of coal-oxygen reaction heat. The temperatures in the three cases at point A are about 360 K, 370 K, and 380 K, respectively. The oxygen concentration is lower than the other two cases. It can be seen that the increase in the coal-oxygen reaction heat accelerates the oxygen consumption rate, thus resulting in a higher spontaneous combustion temperature.

Effect of Activation
Energy. The influences of activation energy E a are investigated in this section. Table 2 lists the activation energy E a values for coal with low, medium, and high cases. 8 Geofluids Figure 9 shows the evolution of oxygen concentration and coal spontaneous combustion temperature at monitoring point A in Cases a, b, and c. Comparing Figure 8 with Figure 9, it can be inferred that the rate of spontaneous combustion of coal is determined by both effects of the rate of coal oxidation and the activation energy E a . The increase in activation energy E a significantly promotes the oxygen consumption rate, while having a sight influence on temperature.

Summary and Conclusions
In this study, we performed a fully coupled dual-porosity and single permeability (DPSP) thermal-hydrological-mechanical-chemical (THMC) modeling of the spontaneous combustion process of underground coal seams, emphasizing the investigation of the influence of pressure difference between oxygen and coal, the rate of coal oxidation, and the activation energy. The main conclusions are drawn based on the simulation results: (1) The simulation results show that the temperature increases with the oxygen propagating into the coal that causes the spontaneous combustion, resulting in the consumption of oxygen along the way and a decrease in gas pressure (2) The expansion of the coal matrix caused by spontaneous combustion of coal reduces the cleat aperture, while the strain induced by matrix shrinkage causes the increase of aperture. The permeability inside the (3) A sensitivity study shows that the increase in pressure difference and coal-oxygen reaction heat contributes to increasing the coal temperature. At the same time, the activation energy has a slight effect on the coal temperature. A more significant pressure difference causes a higher oxygen concentration. A higher value of coal-oxygen reaction heat and activation energy accelerates the oxygen consumption rate, leading to a lower oxygen concentration

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

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