#### Abstract

Combined with the k-*ε* turbulence model of general application, a refined finite element model of a utility tunnel’s gas compartment filled with the methane/air mixture is developed. A series of analyses are made by using the powerful industry-leading computational fluid dynamics (CFD) software flame acceleration simulator (FLACS) to study the shock wave propagation rule in the gas compartment. The longitudinal and transversal distribution laws of the explosion shock wave are gained taking into consideration the spatial characteristics of the gas compartment. The influences of a few parameters, such as initial conditions and section size of the gas compartment, on the shock wave propagation rule are further discussed. The basic procedure for predicting the peak pressure of the blast wave is provided by considering the initial conditions and the gas compartment, and the corresponding injury effect of the explosion wave on the living beings is assessed. The investigation demonstrates that the peak pressure by the coupled effect between the initial conditions is significantly influenced, especially at the upper and lower gas explosion limits. The peak pressure increases gradually as the width or height increases, and both basically meet the linear relation. The proposed method can forecast the peak pressure of the explosion shock wave in the gas compartment accurately. According to the peak pressure longitudinal and transversal distributions of the blast wave, the peak pressure is far greater than the killing pressure threshold in the underground and closed space; consequently, it is not safe for the living beings in the gas compartment.

#### 1. Introduction

Nowadays, the process of the urbanization is accelerated constantly and the scale of the city is also enlarged substantially in the last few decades. The municipal pipeline system becomes more and more complicated and its corresponding management issue becomes increasingly apparent [1, 2]. Accordingly, the utility tunnel, namely, the underground pipe gallery, has emerged in this context [3]. It is an integrated underground pipe ditch or pipeline corridor, which is especially used to lay a variety of municipal pipelines, such as water supply, drainage, heating, gas, electricity, and telecom pipelines. The utility tunnel gradually becomes an important part of the lifeline engineering. It is favored by all walks of life for its easy maintenance, effective use of underground space, good city appearance, and protective effect. The utility tunnel is now extensively applied in major cities around the world [4]. Due to inevitable aging, corrosion, and other human factors, as soon as the leakage accident takes place, it very likely ends in an explosion in the gas compartment in the utility tunnel. Consequently, the research into the propagation law of a gas explosion in the gas compartment has been closely concerned [5, 6].

So far, a range of investigations has been carried out by some experts and scholars from home and abroad to investigate the gas explosion propagation rule in the tunnel. Based on the three-dimensional CFD analysis software AutoReaGas, Pang et al. [7] implemented the numerical simulation analyses to explore the effect of the laneway support spacing on the blast wave properties of a methane/air mixture explosion in a direct laneway and concluded that the support spacing has a great influence on the blast wave distribution. Based on the numerical simulation results from the software AutoReaGas, Zhang et al. [8] put forward a new method to estimate the pressure drop of the explosion wave caused by a premixed methane/air mixture explosion in a closed-ended tunnel and confirm the reliability of this method by the test data. With the aid of the software FLACS, Zhu et al. [9] simulated the explosion process of the methane/air mixture in the direct full-scale tunnel and examined the pressure and its influencing factors, such as the methane gas volume concentration, tunnel blockage ratio, tunnel length, and tunnel cross section. With the application of the large eddy simulation method, Wang et al. [10] investigated the gas explosion shock wave propagation mechanism in the coal mine in consideration of the disaster-causing factors. Combining with the complicated structure characteristics of the goaf, Ke et al. [11] discussed the gas explosion shock wave and flame propagation mechanism in the goaf by using the open-source software OpenFOAM. Wang et al. [12, 13] performed a series of explosion tests to examine the combustion features and explosion characteristics of the methane/ethylene/air mixtures of different equivalence ratios and ethylene volume ratios in a sealed chamber at ambient temperature (20°C) and atmospheric pressure. In addition, they employed a theoretical method based on the adiabatic flame temperature to explore the flammability limits of the methane/air mixtures mixed with the gaseous fuels of different relative humidity. Based on the density functional theory and its detailed mechanism, Su et al. [14] probed intensively the chemical kinetic behavior of the methane/hydrogen mixture thoroughly at the explosion stage and acquired a better grasp and understanding of the methane/hydrogen mixture explosion initiation mechanism. However, the studies regarding the gas explosion propagation rule in the gas compartment are very few. Taking into account the role of the hydrogen, Zhang et al. [15] carried out the finite element analyses to examine the explosion of the methane/hydrogen mixture in the gas compartment of a utility tunnel by means of the commercial CFD analysis software FLACS. From the review of literature mentioned above, the previous works in this field are mostly set in the tunnel of the coal mine and involve its some specific structure characteristics, and lack of the potential risk assessment. In addition, there is an obvious size effect in gas explosion characteristics [16, 17].

Based on the above, and combined with the typical urban gas explosion incidents [18–21] over the latest years, a range of finite element analyses are performed on a typical and representative gas compartment in a utility tunnel to research the propagation rule of the gas blast. A new methodology is presented to forecast the blast wave properties in the gas compartment [22, 23].

#### 2. Gas Compartment

##### 2.1. Brief Introduction of Gas Compartment

A typical and representative gas compartment [24–27] of a real utility tunnel in China’s Beijing City is introduced in this study, as shown in Figure 1, represented by the following parameters: section width, = 2.00 m; section height, *h*_{c} = 4.00 m; length of gas compartment (distance between two adjacent firewalls), *L*_{c} = 200.00 m. The city gas is uniformly mixed together with the air at a specific ratio to create the gas mixture. The gas compartment is filled with the gas mixture at a certain length at one end of it. The gas mixture is ignited by the ignition source, whose temperature is up to about 2000°C. The initial temperature and initial gas volume concentration are 0°C and 9.50%, respectively. The atmospheric pressure is 1.01 × 10^{5} Pa.

The incomplete and complete combustion reaction equations of the gas mixture can be expressed as follows:

##### 2.2. Basic Theory

Computational fluid dynamics is now extensively employed in tackling all sorts of complex issues in fluid mechanics and becomes increasingly significant. It can accurately and reasonably predict a variety of physical phenomena of fluid, such as chemical reaction. After the experimental and theoretical fluid mechanics, the computational fluid dynamic has become a more and more important research tool in the past decades.

###### 2.2.1. Governing Equations

The methane/air mixture explosion in a gas compartment in a utility tunnel can be deemed as a swift and violent combustion process [28, 29]. In the Cartesian (or rectangular) coordinate system, they are provided in partial differential form as follows.

Conservation equation of mass is given bywhere *β*_{ν} is the volume porosity; *β*_{j} is the area porosity in the *j* direction; *ρ* is the mass density; is the mass flow; *V* is the volume; is the velocity in the *j* direction; *x*_{j} is the Cartesian coordinate; and *t* is the time.

Conservation equation of momentum is given bywhere *p* is the hydrostatic pressure; *σ*_{ij} is the stress tensor; *R*_{i} is the distributed resistance in the *x*_{i} direction because of the subgrid obstruction; *R*_{W} is the flow resistance due to the interaction between the fluid and the vessel wall; *ρ*_{0} is the initial mass density; and is the acceleration of gravity in the *i* direction.where *f*_{i} is a constant associated with the obstruction type and the orientation and *A*_{i} = (1 - *β*_{i})/(∆*x*_{i}).

In equation (3), the stress tensor can be given bywhere *μ*_{eff} is the effective viscosity; *k* is the turbulence kinetic energy; and *δ*_{ij} is the Kronecker Delta function.

Enthalpy is an important state parameter characterizing the material system energy in the thermodynamics. The change in the enthalpy is equal to all the energy released by the explosion in the detonation process. Conservation equation of enthalpy:where *h* is the enthalpy; *σ*_{h} is the turbulence constant for the variable ; is the substantial derivative operator; and is the wall heat flux.

Conservation equation of chemical specie mass fraction is given bywhere *m* is the mass of chemical specie; *σ*_{m} is the turbulence constant for the variable *R*_{m}; and *R*_{m} is the reaction rate of the gas.

Additionally, the conservation equations of turbulent kinetic energy and its dissipation rate [30] in the k-*ε* turbulence model can be described as follows:where *k* and *ε* are the turbulent kinetic energy and its dissipation rate, respectively; *σ*_{k} and *σ*_{ε} are the k-*ε* model constants for the variables (*P*_{k} - *β*_{ν}*ρε*) and (C_{1}*P*_{k} - C_{2}*β*_{ν}*ρε*)*ε*/*k*, respectively; *P*_{k} is the turbulent kinetic energy production rate; and *C*_{1} and *C*_{2} are all the constants.

All the conservation equations can be rewritten into a unified compact form and given bywhere *φ* and *ϕ* are the general variables; Γ is the effective turbulent diffusion coefficient; and *S* is the source term.

The expressions for *φ*, *ϕ*, Γ, and *S* are summarized in Table 1.

###### 2.2.2. Other Details

In FLACS, the algorithm to solve the fluid field is the proverbial semi-implicit method for pressure-linked equations (SIMPLE) algorithm in the computational fluid dynamics, which is proposed by Patankar and Spalding [31] in 1972. This algorithm determines the pressure field on the basis of the staggered grid by using the predictor-corrector method and its key issue is to construct the pressure and velocity correction equations.

The classical turbulence models (for instance, *k*-*ε* and RNG turbulence models) are not applicable to the fluid that is located near the wall. The main reason for this is that the viscous force predominates in the fluid. Hence, the introduction of the wall function can improve the modeling of the flow field adjacent to the wall [29]. It is very suitable for the physical field of the gas explosion in the software FLACS.

##### 2.3. Finite Element Model

The grids define the spatial resolution as well as the extents of the simulation volume. The model of the methane/air mixture is illustrated in Figure 2. The grids do not have to be isotropic and equally spaced. Please note that this model is important to have fine enough grids to represent the detailed flow properties and the obstacles properly and the grids must extend quite a large distance from the area of interest to avoid too strong influence from the open boundaries.

The boundary conditions of the numerical simulation domain in Figure 2 are EULER by default. The mass and momentum conservation equations are imposed on the boundary conditions. The boundary conditions should be used with utmost care and attention; otherwise, they can influence the simulated results. The EULER boundary condition may give too high peak pressure in the confined situation.

There is a corresponding stable propagation detonation velocity for a certain gas mixture. The explosion is unstable if the detonation velocity is lower or higher than this stable propagation detonation velocity. The explosion eventually either becomes stable or is over.

#### 3. Model Validation

The gas explosion propagation process in [32, 33] is simulated by using the software FLACS to validate the finite element model.

The parameters representing the square steel pipe in the gas explosion test are as follows: section height and width, *hw*_{p}; length, *l*_{p}; and wall thickness, *t*_{w}. There are two steel pipes made from the 16 Mn steel in this validation study: short pipe (*hw*_{p} = 80 mm; *l*_{p} = 4,000 mm; and = 12 mm) and long pipe (*hw*_{p} = 80 mm; *l*_{p} = 21,000 mm; and = 12 mm). They are filled with the premixed methane/air mixture, and the corresponding gas volume concentration is around 9.50%. The mechanical properties of the pipe material are as follows: yield strength, *f*_{y} = 345 MPa; Young’s modulus, *E* = 2.06 × 10^{11} Pa; Poisson’s ratio, *ν* = 0.30; and mass density, *ρ* = 7,850 kg/m^{3}. The inside of one of them is pasted with a layer of asbestos cloth of 0.80 mm in thickness by using the high temperature resistant adhesive. One end of the steel pipe is closed, while the other end is open. The ignition system is located at the closed end, and the energy of the ignition is around 2 J. The inside of the steel pipe needs to be rubbed down well with the sandpaper so that it can reduce the transition time and distance from deflagration to detonation, thus reducing or even avoiding the effect on the explosion test results. It should be noted that every test is repeated five times.

As listed in Tables 2 and 3, the discrepancies between them are a result of the following four primary reasons. First, the steel pipe is easily deformed slightly owing to the gas blast; thereby it can affect the pressure sensor measure precision. Second, the ignition source is difficult to be well modeled using an ignition temperature [34]. Third, there exists usually a very small amount of the inflammable and explosive gases of other types during the experiment [35]. And fourth, the heat insulation layer formed by the asbestos cloth in the experiment is not the ideal adiabatic boundary condition in the finite element model. But in general, the maximum relative deviations between test and simulation are all less than 5% and are insignificant, and the numerical simulation results coincide very well with the experimental data. Hence, they proved to be accurate and credible in the prediction of the blast wave peak pressures in the short and long pipes.

#### 4. Blast Wave Space Distribution

To evaluate the injury effect on humans and animals and the damage effect on the gas compartment itself, the space distribution [36, 37] of the explosion shock wave in the gas compartment should be examined, which mainly includes the peak pressure longitudinal and transversal distributions. The propagation process of the gas explosion with different filled lengths in the gas compartment is simulated with the application of the reliable commercial CFD software FLACS. The specific research works are as follows.

##### 4.1. Peak Pressure Longitudinal Distribution

On account of the related researches [8], and combined with the numerical results of this study, the blast wave peak pressure along the longitudinal direction can be expressed as follows:where *r*_{0} is the filled length of the gas and *r* is the distance from the measured point to one end of the gas compartment along the longitudinal direction.

The determination coefficients of equations (9) and (10) are 0.9865 and 0.9547, respectively. From equation (10), it can be concluded that the maximum pressure *p*_{Lmax} is 19.75 × (*r*_{0})^{0.84} at *r* = *r*_{0}. At *r*_{0} = 10, 20, 30, 40, and 50 m, they are 138.52, 250.37, 341.83, 438.51, and 525.21 kPa, respectively, as depicted in Figure 3.

As Figure 4 describes, the peak pressure decreases along with the distance away from the gas zone, and the decreasing trend also decelerates. Furthermore, the larger the filled length of the gas is, the higher the peak pressure is. The maximum relative deviations between the simulated and fitting results are all less than 5%, indicating that the fitting and numerical simulation results coincide very well with each other, despite the filled length of the gas.

##### 4.2. Peak Pressure Transversal Distribution

Based on the relevant studies [7], the blast wave peak pressure along the transversal direction can be given as follows:where *zy* is the Cartesian coordinate of the measured point in the *z* or *y* direction and *hw*_{c} is the height or width of the gas compartment.

As Figure 5 illustrates, the peak pressure increases gradually with the increase in |*zy*|; the main reason for this is the reflection [38] of the blast wave that occurs in the internal wall of the gas compartment. Nevertheless, the increasing trend decelerates with the distance away from the gas zone [39]. The peak pressure on the wall of the gas compartment increases by about 22.47% compared with the peak pressure at the center of the section when *r* = *r*_{0}. The maximum relative deviations between the numerical simulation and fitting results are small enough so that the fitting results match well with the simulated results.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

#### 5. Parametric Studies

Unlike the explosive charge explosion, the gas explosion tends to be affected by environmental conditions. Consequently, there is an urgent need to get a comprehensive and systematic understanding of these interfering factors such as initial conditions [40–42] and gas compartment section size [43, 44]. The specific research works are as follows.

##### 5.1. Initial Condition

###### 5.1.1. Initial Temperature

As Figure 6 illustrates, the increase of the initial temperature would increase the gas combustion reaction rate, but decrease the corresponding total mass of the gas, since the gas would be expelled from the gas compartment. It indicates that the peak pressure primarily depends on the total mass of the gas.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

###### 5.1.2. Initial Gas Volume Concentration

As Figure 7 describes, the limited oxygen volume concentration determines the peak pressure. Therefore, as the initial gas volume concentration rises, the peak pressure rises to the maximum when the corresponding initial gas volume is nearly equal to the stoichiometric volume concentration [45] which is around 10%, and then reduces gradually.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

The initial condition influence coefficient can be computed by dividing *p*_{ic} by *p*_{Lmax} from equation (10). Therefore, the difference in the filled length of the gas is eliminated as a matter of course. The corresponding suggested equation can be given by

As Figure 8 presents, the gas explosion parameters are greatly influenced by the initial temperature [46]. The upper gas explosion limit increases but the lower gas explosion limit decreases as the initial temperature increases. Consequently, the peak pressure by the coupled effect between the initial temperature and the initial gas volume concentration is significantly influenced, especially at the upper and lower gas explosion limits [47].

##### 5.2. Gas Compartment Section Size

To examine the effect of the section size on the explosion wave peak pressure, a series of finite element analyses are performed on the gas compartments of different widths and heights, and their results are plotted in Figure 9.

**(a)**

**(b)**

As shown in Figure 9, the peak pressure increases gradually as the width or height increases, and both basically meet the linear relation. Additionally, both width and height enhance the blast wave peak pressure to varying degrees [48], and the height is slightly significant. The main reason for this is the gravity effect.

The section size influence coefficient can be calculated by dividing *p*_{s} by *p*_{Lmax} from equation (10). Therefore, the difference in the filled length of the gas is eliminated for granted. It can be expressed as follows:

An empirical relation among section size influence coefficient, width, and height can be gained using the linear regression analysis and described as equation (13). The determination coefficients of equation (13) are 0.9498. It can be found from Figure 9 that the fitting and finite element analysis results are nearly consistent with each other, despite the width or height of the gas compartment.

#### 6. Suggested Methodology to Estimate the Blast Wave Properties

##### 6.1. Proposed Method Basic Procedure

A new methodology is put forward to forecast the explosion wave peak pressure, whose basic procedures are given as follows:(1)Give gas type, gas filled length, initial condition, gas compartment section size, etc. If not the methane gas, it is converted into the methane equivalent mass based on equation (14) widely employed to calculate the explosion wave peak pressure for different types of explosive charges. The corresponding initial gas volume concentration is also got. Repeat steps (2)–(3). The other gas equivalent mass can be given by where *E*_{Methane} is the methane gas internal energy; *E*_{Other} is the other gas internal energy; and *m*_{Other} is the other gas total mass.(2)The initial condition and section size influence coefficients *γ*_{ic} and *γ*_{s} can be gained on account of equations (12) and (13), respectively. Calculate the maximum pressure *p*_{Lmax} from equation (10). Thus, the maximum pressure considering the initial condition and the gas compartment section size is *γ*_{ic}*γ*_{s}*p*_{Lmax}.(3)The blast wave peak pressure shows an uneven distribution in the gas compartment, and it is an urgent need to obtain the explosion wave peak pressure anywhere. It can be obtained by using equations (9) and (11).(4)Draw the blast wave peak pressure nephogram.

The calculation formulae utilized in the suggested methodology are on account of the former studies and the simulation results. It demonstrates that the proposed methodology can give a good prediction of the blast wave peak pressure in the gas compartment.

##### 6.2. Case Study

To evaluate the explosion shock wave killing power effectively, it is required to develop several quantitative killing criteria, for instance, pressure, impulse, and combination of the above two. Reference [23] provides 45.7 kPa for eardrum rupture with a failure probability of about 10% and 103.4 kPa in the duration of lung damage by 50 ms. In the meantime, the analogical killing criteria are to be found in other countries. In the research background of a gas compartment of a utility tunnel in China’s Beijing City, the propagation process of the gas explosion is simulated by using the powerful industry-leading CFD software FLACS to analyze and discuss the blast wave space distribution.

As Figure 10 illustrates, the longer the filled length of the gas is and the lower the initial temperature is, the higher the blast wave peak pressure becomes. As the initial gas volume concentration increases, the peak pressure increases firstly and decreases afterward slowly.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(g)**

**(h)**

**(i)**

According to the longitudinal and transversal distributions, the peak pressure is far greater than the killing pressure threshold in the underground and closed space; consequently, it is not safe for the living beings in the gas compartment.

#### 7. Conclusion

In this work, a new methodology has been put forward to forecast the gas explosion shock wave peak pressure in the gas compartment. A range of analyses has been carried out on this issue. Some conclusions are deduced as follows:(1)The peak pressure decreases along with the distance away from the gas zone, and the decreasing trend also decelerates; moreover, the longer the filled length of the gas is, the higher the peak pressure is.(2)The peak pressure by the coupled effect between the initial temperature and the initial gas volume concentration is significantly influenced, especially at the upper and lower gas explosion limits. The peak pressure increases gradually as the width or height increases, and both basically meet the linear relation.(3)The analysis procedure which is complete and rigorous demonstrates that the proposed methodology can effectively estimate the explosion wave peak pressure in the gas compartment.(4)According to the longitudinal and transversal distributions, the peak pressure is far greater than the killing pressure threshold in the underground and closed space; consequently, it is not safe for the living beings in the gas compartment.

#### Data Availability

The simulation data used to support the findings of this study are included within the article.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.

#### Acknowledgments

The authors gratefully acknowledge the financial support provided by the Key Research and Development (R&D) Program of Tangshan (no. 19150232E), Xijing University's Innovation and Entrepreneurship Training Program for College Students in 2021 (X202112715007), and Postdoctoral Research Project of Chongqing (no. Xm2017189).