Using CFD to Simulate the Concentration of Polluting and Harmful Gases in the Roadway of Non-Metallic Mines Reveals Its Migration Law

: The green and pollution-free mining of resources has always been a research ﬁeld that people have focused on. In the process of mining resources, the production of CO, SO 2 and other pollutants directly affects the health of miners and the atmospheric environment in the mining area. Therefore, it is particularly important to deal with and control the polluting gases generated by mining. Taking the underground roadway of a coal mine in Hengdong City, Hunan Province, as the research object, we studied the migration law of pollutant gas there. Comsol software was used to determine the changing state of pollutant gas migration in the roadway, and a simulation model of the wind ﬁeld and the pollutant concentration ﬁeld in the roadway under turbulent conditions was established. The results showed that, when the air ﬂow moved to the front face of the roadway, it generated backﬂow to form a counterclockwise-rotating air ﬂow vortex. Here, the air ﬂow stagnated, hindering the diffusion of pollutants. The gas moved with the air ﬂow in the roadway, and the ﬂow’s velocity decreased in the middle of the roadway, causing pollutants to accumulate. The whole wind ﬁeld tended to be stable at a plane 25 m from the roadway’s outlet. This indicates that the middle part of the roadway is the place where the polluted gas accumulates, and it is of representative signiﬁcance to study the concentration of the polluted gas in the roadway in this section. signiﬁcantly in the ﬁrst 3 min. The concentration of polluting gas in the upper part of this section was clearly higher than that in other parts at 0–20 m. The concentration of polluting gas in the upper part of this section was still the highest at 20–40 m, and the difference was smaller than that at 0–20 m. As the time increased, the overall concentration of pollutant gas decreased and the difference in the local concentration decreased.


Introduction
In recent years, the development of the economy and society has led to a continuous improvement in the degree of mechanization of mining and a continuous increase in mining intensity [1,2]. The open-pit mining method, among the common mining methods, destroys the surface soil and the ecological balance of animals and plants, and the damage to the surrounding natural environment is direct and obvious. By comparison, underground mining of ores causes less damage to the surrounding natural environment [3][4][5][6]. However, the pollution caused by ventilation in underground mines is still unavoidable [2]. Dust and fugitive dust generated by the operation of large-scale machinery and blasting operations generate a large number of gaseous pollutants such as carbon monoxide, nitrogen oxides and sulfur dioxide, and a range of pollution problems. When the pollutant concentration reaches a certain level, it directly damages the normal physiological function of underground workers [7][8][9][10]. The pollutants generated from underground operations are physically transported along the mine's ventilation network to reduce the harmful concentrations. At the same time, the closure and excavation of some working roadways force the ventilation network's structure to change, causing ventilation defects, which make the accumulation of local polluting gases reach harmful concentrations, endangering the health of operators [11][12][13][14]. For example, the Doyon Westwood mine in Abitiby, Quebec, has sulfide pollution problems, while the Pulang copper mine in Yunnan has carbon monoxide pollution problems [14,15]. At present, the research on toxic gases in coal mines mainly focuses on the amount of indicator gas generated by the spontaneous combustion of coal. Through analysis of the change trend in the concentration of several special gas products generated by the spontaneous combustion and oxidation of coal, it is possible to determine whether there is likely to be a spontaneous combustion disaster in coal mines [16,17]. Metal mines generally pay attention to the migration laws and the treatment and disposal of the dust and toxic and harmful gases generated by blasting. The temporal and spatial evolution of the characteristics of CO gas in blasting smoke have been studied and analyzed under specific working conditions in order to reduce the number of pollutants at the source of blasting generated by different blasting methods [18,19]. Some non-uranium mines such as lead-zinc ore, tungsten ore and other non-ferrous metal ores are often contaminated by radioactive uranium, thorium and other nuclides, which will decay during the mining process, resulting in pollution by radioactive radon gas and occupational diseases such as lung cancer [20][21][22].
The internal structure of an underground mine is complex. In order to provide a suitable environment for workers to work, it is necessary to lay a ventilation pipe network to strengthen the capacity for ventilation [23]. However, the polluting gas generated by the blasting operation for tunnel excavation moves in the tunnel with the wind flow, which not only pollutes the surrounding environment but also endangers the health of workers [24]. Wang et al. studied the dust pollution generated by blasting and construction mining in open-pit coal mines, and summed up the corresponding accident prevention measures and optimal strategies for mine design [25]. Harris analyzed the CO generated by mine-blasting operations and showed that applying negative pressure to the borehole was the most effective method of removing CO from the surface and minimizing its accumulation and migration [26]. Rudakov et al. improved the filling of explosives and provided a method to reduce nitrogen oxide emissions by using highly active catalysts as a part of the molding and tamping process [27]. The results showed that using zinc carbonate (ZnCO 3 ) as a filler reduced nitrogen oxide emissions by 40% without the addition of catalysts. Adhikari et al. studied the composition and retention characteristics of blasting flue gas, and proposed effective measures to ensure workers' health as well as clean production in wells [28].
With the replacement of various CPUs and the rapid development of computer technology, the methods of computational fluid dynamics theory have been applied in various numerical simulation software packages, and the use of numerical simulation software to study fluid laws has become the mainstream today [29]. Feng et al. used Fluent software to simulate the diffusion behavior of harmful gases after tunnel blasting in a high-altitude environment, and revealed the dynamic law that the equivalent mass concentration of CO increases exponentially with an increase in altitude [30]. Xie et al. used the computational fluid dynamics (CFD) method to analyze the characteristics of the tunnel flow field and established a fluid model after blasting, as well as determining the area where personnel were more seriously exposed to dust [31]. Torno used the CFD method to establish a 4D dynamic model of flue gas, which accurately revealed the gas's dilution behavior [32]. Sa et al. used Fluent software to simulate the change in the concentration of ventilation dust in the cavern stopes and the dust trajectory under different wind speeds within 20 min after blasting. The distribution law of the dust concentration in cavern stopes and the influence of the wind speed on dust were accurately revealed [33]. Huang et al. used the CFD method to study the characteristics of the temporal and spatial evolution of CO under different working conditions, and applied the Grey relational analysis method to study the correlation values of factors affecting the diffusion coefficient of CO in the roadway [34]. Xie et al. used CFD analysis to evaluate the optimal angle of air supply, which has practical significance for improving dust removal and reducing the energy consumption of ventilation [35]. Wang et al. used the CFD method to establish the evolution of the process of explosion-dust diffusing, which provided a theoretical basis for rationally arranging the production of deep-pit open-cast mines and avoiding air pollution [36]. Seung-Chul Lee et al. used CFD analysis to evaluate the ventilation characteristics of fans with different Sustainability 2022, 14, 13349 3 of 16 layout elevation angles and determined the optimal angle, which has practical significance for accelerating the diffusion rate of polluting gases [37]. Roberto Bubbico et al. used the CFD method to study the evolution process of gas caused by the release of toxic substances and established a turbulent flow model for the diffusion of toxic substances within 5 min of the start of the release of the toxic substances [38]. Ping Chang et al. used the CFD method to analyze the airflow characteristics of the working face and the distribution of DPM concentrations, and established the corresponding fluid model [39]. Menéndez et al. used a one-dimensional mathematical model and a three-dimensional CFD numerical simulation to determine the safe re-entry time after tunnel blasting, and used gas sensors to measure the consistency of the CFD and the three-dimensional numerical simulations [40].
At present, scholars have carried out much research on gas pollution in underground mine tunnels, but research on the concentration and migration law of polluting and harmful gases in specific nonmetallic mine tunnels needs to be deepened further. It is still effective and reliable to use numerical simulation software to study the concentration and migration law of polluting and harmful gases in roadways. Therefore, exploration of the migration law of gas concentration in the roadway can provide a scientific basis for the formulation of measures for the prevention and control of pollutant gas in the later stages and to ensure the health of operators [41].

Construction of a Numerical Simulation Model of Polluting Gas in the Roadway
As a branch of fluid mechanics, CFD is widely used in studies into the transport of polluting and harmful gases. The toxic gases generated by mining operations flow in the tunnels and along the air passages in underground mines, causing the accumulation of polluting gases in abnormal ventilation areas. Because of the irregular distribution of underground roadways, the supporting ventilation system and the monitoring system are greatly affected by external conditions, and the gas accumulation phenomenon caused by them is difficult to detect initially. CFD can analyze and predict the pollution and harmful gases in the roadway through numerical simulations. The change in the flow field can accurately reveal the transport law of the polluting gas in the roadway [42]. COMSOL Multiphysics software was used to numerically simulate the polluting gas's flow field in the actual roadway, to quickly and effectively obtain the change in the flow field's status, to determine the area of gas accumulation, and provide scientific guidance for the subsequent prevention of pollution gas and for environmental protection [43].

Concentration Diffusion Control Equation
The fluid flow in an underground mine tunnel obeys the Navier-Stokes equation, and must comply with the general laws of mechanical motion and the fluid mechanics laws of mass conservation, momentum conservation and energy conservation. The equation governing the concentration and diffusion is as follows [44].
(1) Continuity equation where ρ is the gas density, and Vi represents the velocity (Vx, Vy, Vz) in the directions of the coordinate axes x, y and z. A fluid follows the law of conservation of mass during the process of flowing [45]. A representation of the gas micro-element in the Cartesian coordinate system is shown in Figure 1. (2) Momentum conservation equation The various forces exerted on the micro-element by the outside world act together on the movement of the micro-element, so that the micro-element changes with time. The equation for its conservation of momentum is expressed as: where is represents the fluid's dynamic viscosity, G represents gravitational acceleration, P represents the absolute pressure and P  represents the density of air.
(3) Energy conservation equation The energy transfer in the process of the flow of polluting gas in the roadway follows the first law of thermodynamics, and its expression is as follows [46]: where Cp is the specific heat capacity, in J/(kg·°C); T is the temperature; k is the fluid heat transfer coefficient, in W/(m 2 ·°C); and ST is the viscous dissipation term.
(4) Component equation The diffusion of gas in the atmospheric environment needs to satisfy the law of conservation of the components' mass, and the expression is as follows where CZ is the concentration volume of the component Z; where is the mass concentration of the component Z, in g/m 3 ; DZ is the diffusion coefficient of the component Z and Fz is the production rate of the chemical reaction, in g/(s·m 3 ).
In addition, the Reynolds-averaged Navier-Stokes method was used to study the variation law of the pollutant gas's flow field in the tunnel.

Model Site Selection and Construction of the Geometric Model
The mine under study is located near Xiaoshanchong, northeast of Hengdong County, Hengyang City, Hunan Province. The geographic coordinates have a longitude range of 113.073164° to 113.080134° and a latitude range of 27.139523° to 27.137665°. The location is shown in Figure 2. A certain roadway was selected as the research object. This roadway is a single-ended roadway. Because of the roadway's structural characteristics,  (2) Momentum conservation equation The various forces exerted on the micro-element by the outside world act together on the movement of the micro-element, so that the micro-element changes with time. The equation for its conservation of momentum is expressed as: where is represents the fluid's dynamic viscosity, G represents gravitational acceleration, P represents the absolute pressure and P α represents the density of air.
(3) Energy conservation equation The energy transfer in the process of the flow of polluting gas in the roadway follows the first law of thermodynamics, and its expression is as follows [46]: where C p is the specific heat capacity, in J/(kg· • C); T is the temperature; k is the fluid heat transfer coefficient, in W/(m 2 · • C); and S T is the viscous dissipation term.
(4) Component equation The diffusion of gas in the atmospheric environment needs to satisfy the law of conservation of the components' mass, and the expression is as follows where C Z is the concentration volume of the component Z; where is the mass concentration of the component Z, in g/m 3 ; DZ is the diffusion coefficient of the component Z and F z is the production rate of the chemical reaction, in g/(s·m 3 ). In addition, the Reynolds-averaged Navier-Stokes method was used to study the variation law of the pollutant gas's flow field in the tunnel.

Model Site Selection and Construction of the Geometric Model
The mine under study is located near Xiaoshanchong, northeast of Hengdong County, Hengyang City, Hunan Province. The geographic coordinates have a longitude range of 113.073164 • to 113.080134 • and a latitude range of 27.139523 • to 27.137665 • . The location is shown in Figure 2. A certain roadway was selected as the research object. This roadway is a single-ended roadway. Because of the roadway's structural characteristics, the airflow in the roadway is often unsmooth. When the concentration of CO, SO 2 and other toxic polluting gases in the road is too high, additional means are needed to assist the diffusion of the polluting toxic gas in the roadway. Therefore, diffusion of the concentrations of polluting and harmful gas occurs in the exhaust air of the tunnel. In the area where the gas flow rate is abnormal and the gas accumulates in the tunnel, there is also an excessively high concentration of polluting gases, which has a significant impact on the health of the operators and the tunnel environment. the airflow in the roadway is often unsmooth. When the concentration of CO, SO2 and other toxic polluting gases in the road is too high, additional means are needed to assist the diffusion of the polluting toxic gas in the roadway. Therefore, diffusion of the concentrations of polluting and harmful gas occurs in the exhaust air of the tunnel. In the area where the gas flow rate is abnormal and the gas accumulates in the tunnel, there is also an excessively high concentration of polluting gases, which has a significant impact on the health of the operators and the tunnel environment. The selected roadway is 40 m below the ground surface. Roadway L09 was selected as the simulation object; the size of its cross-section is 4 m × 3.2 m, its area is about 11.2 m 2 and the total length of the roadway is 45 m. This section has the characteristics of a narrow roadway, with a long excavation length and a low ventilation pressure. In order to better simulate the flow law of the polluting gas generated by blasting, the assumptions of the simulation were as follows. (1) The air in the roadway before the blasting operation was ignored. (2) The gas generated by the blasting was a three-dimensional incompressible gas. (3) The temperature in the roadway was constant. (4) The influence of other equipment on ventilation was ignored. (5) The flexible fireproof material of the ventilation equipment was simplified, the equipment's radius was 0.4 m, the inlet's wind speed was 10 m/s and the distance from the tunnel face was 12 m. (6) The modeled outlet was set as a free outlet boundary condition, where the outlet's pressure was 1 atm, the roughness constant was 0.7, and the roughness of particle height was 0.3. (7) The fan blew air to this roadway from the fresh air flow, bringing out the residual polluting poisonous gas.
Based on these basic conditions, a geometric model (shown in Figure 3) was established. The fluid properties were determined according to the actual field conditions and theoretical calculations. The selected roadway is 40 m below the ground surface. Roadway L09 was selected as the simulation object; the size of its cross-section is 4 m × 3.2 m, its area is about 11.2 m 2 and the total length of the roadway is 45 m. This section has the characteristics of a narrow roadway, with a long excavation length and a low ventilation pressure. In order to better simulate the flow law of the polluting gas generated by blasting, the assumptions of the simulation were as follows. (1) The air in the roadway before the blasting operation was ignored. (2) The gas generated by the blasting was a three-dimensional incompressible gas.
(3) The temperature in the roadway was constant. (4) The influence of other equipment on ventilation was ignored. (5) The flexible fireproof material of the ventilation equipment was simplified, the equipment's radius was 0.4 m, the inlet's wind speed was 10 m/s and the distance from the tunnel face was 12 m. (6) The modeled outlet was set as a free outlet boundary condition, where the outlet's pressure was 1 atm, the roughness constant was 0.7, and the roughness of particle height was 0.3. (7) The fan blew air to this roadway from the fresh air flow, bringing out the residual polluting poisonous gas.
Based on these basic conditions, a geometric model (shown in Figure 3) was established. The fluid properties were determined according to the actual field conditions and theoretical calculations.

Grid Generation and Verification of Iterative Convergence
The geometric model was meshed by using construction mesh technology. The total number of elements was 140,999, and the minimum element quality was 0.05924. The grid was surrounded by fixed boundaries, and the entrance and exit were free-exit boundary conditions. The results of the mesh division are shown in Figure 4.

Grid Generation and Verification of Iterative Convergence
The geometric model was meshed by using construction mesh technology. The total number of elements was 140,999, and the minimum element quality was 0.05924. The grid was surrounded by fixed boundaries, and the entrance and exit were free-exit boundary conditions. The results of the mesh division are shown in Figure 4.  By setting the wind speed at the inlet of the ventilation pipe in the model to 10 m/s, the pressure outlet conditions, the boundary conditions and the movement of polluting gas in the roadway were analyzed. After setting the relevant physical parameters of the flow field, CFD software was used to start the iterative calculation. A diagram of the iter-

Grid Generation and Verification of Iterative Convergence
The geometric model was meshed by using construction mesh technology. The total number of elements was 140,999, and the minimum element quality was 0.05924. The grid was surrounded by fixed boundaries, and the entrance and exit were free-exit boundary conditions. The results of the mesh division are shown in Figure 4.  By setting the wind speed at the inlet of the ventilation pipe in the model to 10 m/s, the pressure outlet conditions, the boundary conditions and the movement of polluting gas in the roadway were analyzed. After setting the relevant physical parameters of the flow field, CFD software was used to start the iterative calculation. A diagram of the iter- By setting the wind speed at the inlet of the ventilation pipe in the model to 10 m/s, the pressure outlet conditions, the boundary conditions and the movement of polluting gas in the roadway were analyzed. After setting the relevant physical parameters of the flow field, CFD software was used to start the iterative calculation. A diagram of the iterative convergence is shown in Figure 5. This study compared and analyzed the four convergence standards of 10 0 , 10 −1 , 10 −2 , 10 −3 and 10 −4 . It can be seen that the convergence standards of 10 0 to 10 −2 differed greatly and did not represent convergence. The difference between the results of the 10 −3 and 10 −4 convergence standards was small;therefore, the simulated convergence criterion was 10 −4 . Each time step in the diagram awas in a convergence state. The solver reached the given convergence standard after 85 iterations, which indicates that the numerical simulation was effective. vergence standards of 10 0 , 10 −1 , 10 −2 , 10 −3 and 10 −4 . It can be seen that the convergence standards of 10 0 to 10 −2 differed greatly and did not represent convergence. The difference between the results of the 10 −3 and 10 −4 convergence standards was small;therefore, the simulated convergence criterion was 10 −4 . Each time step in the diagram awas in a convergence state. The solver reached the given convergence standard after 85 iterations, which indicates that the numerical simulation was effective. The actual wind speed and temperature data were monitored by three detectors arranged at 30 m and 15 m along the X axis of the roadway at an interval of 1 min.

Setting of the Detection Points
According to the actual physical parameters of the roadway, two detection sections were selected as the data collection, which were located in the middle and last two sections of the roadway, with a distance of 15 m and a total of 6 detection points, surfaced. We used CAD software to generate the positions in the roadway corresponding to each detection point. The first detection section was 15 m away from the tunnel face, as shown in Figure 6. The test results are shown in Table 1.  The actual wind speed and temperature data were monitored by three detectors arranged at 30 m and 15 m along the X axis of the roadway at an interval of 1 min.

Setting of the Detection Points
According to the actual physical parameters of the roadway, two detection sections were selected as the data collection, which were located in the middle and last two sections of the roadway, with a distance of 15 m and a total of 6 detection points, surfaced. We used CAD software to generate the positions in the roadway corresponding to each detection point. The first detection section was 15 m away from the tunnel face, as shown in Figure 6. The test results are shown in Table 1.
between the results of the 10 −3 and 10 −4 convergence standards was small;therefore, the simulated convergence criterion was 10 −4 . Each time step in the diagram awas in a convergence state. The solver reached the given convergence standard after 85 iterations, which indicates that the numerical simulation was effective. The actual wind speed and temperature data were monitored by three detectors arranged at 30 m and 15 m along the X axis of the roadway at an interval of 1 min.

Setting of the Detection Points
According to the actual physical parameters of the roadway, two detection sections were selected as the data collection, which were located in the middle and last two sections of the roadway, with a distance of 15 m and a total of 6 detection points, surfaced. We used CAD software to generate the positions in the roadway corresponding to each detection point. The first detection section was 15 m away from the tunnel face, as shown in Figure 6. The test results are shown in Table 1.

Data Collection and Processing
Due to the complex construction environment on site, the detection data were affected by environmental restrictions and human factors. The internal storage function of the portable gas detector was used to save the detection data, import them into the computer through the data cable and then clean the data in the dataset. The outlier detection method based on the local outlier factor (LOF) algorithm was used to detect the abnormal values of the acquired historical downhole gas concentrations. If an abnormal value was detected, the average value of the three concentration data points before the abnormal value was used to replace the abnormal value, as shown in Figure 7. The local outlier factor of the data (P) was determined by calculating the local reachable density of the P for the historical gas concentration data obtained by the LOF algorithm. When the local outlier factor was less than or equal to 1, the point P was not an abnormal point in the neighborhood or even a reunion point. When this ratio was greater than 1, the distance between point P and the points in its neighborhood was large, and P may be an abnormal point. When calculating the local reachable density of data point P, we first determined the k-th distance of P, where the k-th distance was the distance from the data point farthest from P to P in P's k-neighborhood. The k-neighborhood of point P was determined by the following definitions:  There are at least k points O' excluding P in the set, with the result that the distance from O' to P is less than or equal to the kth distance of P.

Data Collection and Processing
Due to the complex construction environment on site, the detection data were affected by environmental restrictions and human factors. The internal storage function of the portable gas detector was used to save the detection data, import them into the computer through the data cable and then clean the data in the dataset. The outlier detection method based on the local outlier factor (LOF) algorithm was used to detect the abnormal values of the acquired historical downhole gas concentrations. If an abnormal value was detected, the average value of the three concentration data points before the abnormal value was used to replace the abnormal value, as shown in Figure 7.

Data Collection and Processing
Due to the complex construction environment on site, the detection data were affected by environmental restrictions and human factors. The internal storage function of the portable gas detector was used to save the detection data, import them into the computer through the data cable and then clean the data in the dataset. The outlier detection method based on the local outlier factor (LOF) algorithm was used to detect the abnormal values of the acquired historical downhole gas concentrations. If an abnormal value was detected, the average value of the three concentration data points before the abnormal value was used to replace the abnormal value, as shown in Figure 7. The local outlier factor of the data (P) was determined by calculating the local reachable density of the P for the historical gas concentration data obtained by the LOF algorithm. When the local outlier factor was less than or equal to 1, the point P was not an abnormal point in the neighborhood or even a reunion point. When this ratio was greater than 1, the distance between point P and the points in its neighborhood was large, and P may be an abnormal point. When calculating the local reachable density of data point P, we first determined the k-th distance of P, where the k-th distance was the distance from the data point farthest from P to P in P's k-neighborhood. The k-neighborhood of point P was determined by the following definitions: The local outlier factor of the data (P) was determined by calculating the local reachable density of the P for the historical gas concentration data obtained by the LOF algorithm. When the local outlier factor was less than or equal to 1, the point P was not an abnormal point in the neighborhood or even a reunion point. When this ratio was greater than 1, the distance between point P and the points in its neighborhood was large, and P may be an abnormal point. When calculating the local reachable density of data point P, we first determined the k-th distance of P, where the k-th distance was the distance from the data point farthest from P to P in P's k-neighborhood. The k-neighborhood of point P was determined by the following definitions:

•
There are at least k points O' excluding P in the set, with the result that the distance from O' to P is less than or equal to the kth distance of P. • There are at most k − 1 points O' in the set excluding the point P, with the result that the distance from O' to P is smaller than the kth distance of P.
For example, in the example shown in Figure 7, k = 2 for the point X(1), and the two closest points to it are the point X(1) and the point X(4); in other words, the 2-neighborhood of X(1) includes point X(1) and point X(4).
We then determined the reachable distance between each data point O and P in the k-neighborhood of P, where the reachable distance between O and P was the maximum of the distance between O and P and the kth distance of P. After that, the local reachability density of P could be determined according to the following Formula (5): where Ird k (P) represents the local reachability density of P, N k (P) represents the k-neighborhood of P and the reachable distance between O and P. If P and the points in the neighborhood are of the same class, the reachable distance is more likely to be smaller, resulting in a larger local reachability density. If both P and the points in the neighborhood are far away, the reachable distance may be large, and the local reachable density will be small. When the local outlier factor of the data is determined according to the local reachability density, the local outlier factor of the data can be calculated according to the following Formula (6):

Variation Law of the Wind Field under Specific Conditions
We selected the three-direction section of the diagram of the results of the simulation model, in which the Z-axis direction was selected as a plane with a height of 2.2 m, and the Y-axis direction was selected as a plane with a width of 2 m, as shown in Figure 8. We selected planes with lengths of 10 m, 20 m, 30 m, 35 m and 40 m in the X-axis direction to analyze the distribution of the wind flow field in the roadway, as shown in Figure 9. • There are at least k points O' excluding P in the set, with the result that the distance from O' to P is less than or equal to the kth distance of P.

•
There are at most k − 1 points O' in the set excluding the point P, with the result that the distance from O' to P is smaller than the kth distance of P.
For example, in the example shown in Figure 7, k = 2 for the point X(1), and the two closest points to it are the point X(1) and the point X(4); in other words, the 2-neighborhood of X(1) includes point X(1) and point X(4).
We then determined the reachable distance between each data point O and P in the k-neighborhood of P, where the reachable distance between O and P was the maximum of the distance between O and P and the kth distance of P. After that, the local reachability density of P could be determined according to the following Formula (5): where Irdk(P) represents the local reachability density of P, Nk(P) represents the k-neighborhood of P and the reachable distance between O and P. If P and the points in the neighborhood are of the same class, the reachable distance is more likely to be smaller, resulting in a larger local reachability density. If both P and the points in the neighborhood are far away, the reachable distance may be large, and the local reachable density will be small. When the local outlier factor of the data is determined according to the local reachability density, the local outlier factor of the data can be calculated according to the following Formula (6):

Variation Law of the Wind Field under Specific Conditions
We selected the three-direction section of the diagram of the results of the simulation model, in which the Z-axis direction was selected as a plane with a height of 2.2 m, and the Y-axis direction was selected as a plane with a width of 2 m, as shown in Figure 8. We selected planes with lengths of 10 m, 20 m, 30 m, 35 m and 40 m in the X-axis direction to analyze the distribution of the wind flow field in the roadway, as shown in Figure 9.  It can be seen from the cross-sections of the air flow field inside the roadway in Figures 8 and 9 that the air duct transmitted the flow of fresh air to the front section of the roadway at a speed of 10 m/s, and the wind speed was the highest at the exit. When the length was 10 to 30 m, the wind speed at the lower right-hand side of the roadway was high. The airflow moved forward along the roadway wall, the flow distance of the affected area increased and expanded continuously, and the wind speed also decreased gradually. When the airflow moved to the front face of the roadway, it quickly generated a backflow phenomenon, forming an anticlockwise airflow vortex. The wind flow in the vortex was stagnant. When the polluting gas molecules were in the range of the vortex, they rotated with the wind field, which was not conducive to the diffusion of pollutants, and the polluting gas accumulated easily. It can be seen that the wind speed was the greatest within 5-15 m of the front section of the roadway and gradually decreased after that. The whole wind field tended to be stable in the section 25 m away from the roadway's exit. The measured wind speed in Table 1 was close to the results of the simulation, which proves that the simulation's results are reliable and can be used as a conclusion.

Change Law of Polluting Gas under the Action of the Wind Field
After the simulation calculation was completed and the conclusions drawn from the simulation were verified to be reliable, the migration law of the polluting gas in the roadway under the action of the wind field was studied. The carbon monoxide gas component was added to the simulation for transient research. Lines 1, 5, 10 and 15 are shown in Figure 10. The concentration field was the lowest at 15 min. It can be seen that the pollutant gas was reduced by the air flow of the roadway, and the pollutant gas was essentially discharged 15 min later. Because of the eddy current in the wind field, the pollutants in the left-and right-hand corners of the roadway were not easily spread, and the flow of the pollutant gas in the 25 m plane tended to be stable. It can be seen from the cross-sections of the air flow field inside the roadway in Figures 8 and 9 that the air duct transmitted the flow of fresh air to the front section of the roadway at a speed of 10 m/s, and the wind speed was the highest at the exit. When the length was 10 to 30 m, the wind speed at the lower right-hand side of the roadway was high. The airflow moved forward along the roadway wall, the flow distance of the affected area increased and expanded continuously, and the wind speed also decreased gradually. When the airflow moved to the front face of the roadway, it quickly generated a backflow phenomenon, forming an anticlockwise airflow vortex. The wind flow in the vortex was stagnant. When the polluting gas molecules were in the range of the vortex, they rotated with the wind field, which was not conducive to the diffusion of pollutants, and the polluting gas accumulated easily. It can be seen that the wind speed was the greatest within 5-15 m of the front section of the roadway and gradually decreased after that. The whole wind field tended to be stable in the section 25 m away from the roadway's exit. The measured wind speed in Table 1 was close to the results of the simulation, which proves that the simulation's results are reliable and can be used as a conclusion.

Change Law of Polluting Gas under the Action of the Wind Field
After the simulation calculation was completed and the conclusions drawn from the simulation were verified to be reliable, the migration law of the polluting gas in the roadway under the action of the wind field was studied. The carbon monoxide gas component was added to the simulation for transient research. Lines 1, 5, 10 and 15 are shown in Figure 10. The concentration field was the lowest at 15 min. It can be seen that the pollutant gas was reduced by the air flow of the roadway, and the pollutant gas was essentially discharged 15 min later. Because of the eddy current in the wind field, the pollutants in the left-and right-hand corners of the roadway were not easily spread, and the flow of the pollutant gas in the 25 m plane tended to be stable. In order to analyze the change in the concentration of the polluting gas in the section, it was necessary to study the distribution of the concentration of each section of the roadway. The distribution of concentration at 1, 3, 6, 9, 12, and 15 min in a plane with a length of 40 m was selected, as shown in Figure 11. The overall concentration of polluting gas in this section was low at 12 to 15 min, and the concentration of polluting gas increased significantly after 15 min. The concentration of pollutant gas at the top of the section was the highest, followed by the two sides. The concentration of the pollutant gas at the top and the two sides of the section was always higher than that in the middle, and the concentration of pollutant gas in the middle was the lowest.  In order to analyze the change in the concentration of the polluting gas in the section, it was necessary to study the distribution of the concentration of each section of the roadway. The distribution of concentration at 1, 3, 6, 9, 12, and 15 min in a plane with a length of 40 m was selected, as shown in Figure 11. The overall concentration of polluting gas in this section was low at 12 to 15 min, and the concentration of polluting gas increased significantly after 15 min. The concentration of pollutant gas at the top of the section was the highest, followed by the two sides. The concentration of the pollutant gas at the top and the two sides of the section was always higher than that in the middle, and the concentration of pollutant gas in the middle was the lowest. In order to analyze the change in the concentration of the polluting gas in the section, it was necessary to study the distribution of the concentration of each section of the roadway. The distribution of concentration at 1, 3, 6, 9, 12, and 15 min in a plane with a length of 40 m was selected, as shown in Figure 11. The overall concentration of polluting gas in this section was low at 12 to 15 min, and the concentration of polluting gas increased significantly after 15 min. The concentration of pollutant gas at the top of the section was the highest, followed by the two sides. The concentration of the pollutant gas at the top and the two sides of the section was always higher than that in the middle, and the concentration of pollutant gas in the middle was the lowest.  The distribution of concentration at 0, 1, 3, 6, 9, 12 and 15 min in the plane with a width of 2 m is shown in Figure 12. The concentration of polluting gas in the upper part of this section was the highest compared with other parts, and the concentration changed significantly in the first 3 min. The concentration of polluting gas in the upper part of this section was clearly higher than that in other parts at 0-20 m. The concentration of polluting gas in the upper part of this section was still the highest at 20-40 m, and the difference was smaller than that at 0-20 m. As the time increased, the overall concentration of pollutant gas decreased and the difference in the local concentration decreased. The distribution of concentration at 0, 1, 3, 6, 9, 12 and 15 min in the plane with a width of 2 m is shown in Figure 12. The concentration of polluting gas in the upper part of this section was the highest compared with other parts, and the concentration changed significantly in the first 3 min. The concentration of polluting gas in the upper part of this section was clearly higher than that in other parts at 0-20 m. The concentration of polluting gas in the upper part of this section was still the highest at 20-40 m, and the difference was smaller than that at 0-20 m. As the time increased, the overall concentration of pollutant gas decreased and the difference in the local concentration decreased.  The distribution of concentration at 0, 1, 3, 6, 9, 12 and 15 min in the plane with a height of 2 m is shown in Figure 13. The concentration of pollutant gas generally decreased. The concentration of polluting gas at the inlet and upper part of this section was relatively high. Similarly, the concentration changed significantly in the first three minutes. The concentration field at 25 m tended to be stable. The concentration of polluting gas decreased with time. After 15 min, the concentration could be ignored, and the polluting gas had largely emptied. The distribution of concentration at 0, 1, 3, 6, 9, 12 and 15 min in the plane with a height of 2 m is shown in Figure 13. The concentration of pollutant gas generally decreased. The concentration of polluting gas at the inlet and upper part of this section was relatively high. Similarly, the concentration changed significantly in the first three minutes. The concentration field at 25 m tended to be stable. The concentration of polluting gas decreased with time. After 15 min, the concentration could be ignored, and the polluting gas had largely emptied.
On the basis of the distribution of concentration shown in Figures 10-13, it can be seen that the concentration of carbon monoxide moved with the ventilation, and the highest concentration of polluting gas appeared in the upper part of the roadway. Turbulence in the concentration field at the length of 40 m made it difficult for pollutants in the left-and right-hand corners of the roadway to diffuse and ventilate. The concentration field at a length of 25 m tended to be stable. This plane was selected as the detection surface for the concentration of polluting gases, and the concentration values of polluting gases on this plane were collected to provide an effective scientific basis for predicting the concentration of polluting gases in subsequent sections.  On the basis of the distribution of concentration shown in Figures 10-13, it can be seen that the concentration of carbon monoxide moved with the ventilation, and the highest concentration of polluting gas appeared in the upper part of the roadway. Turbulence in the concentration field at the length of 40 m made it difficult for pollutants in the leftand right-hand corners of the roadway to diffuse and ventilate. The concentration field at a length of 25 m tended to be stable. This plane was selected as the detection surface for the concentration of polluting gases, and the concentration values of polluting gases on this plane were collected to provide an effective scientific basis for predicting the concentration of polluting gases in subsequent sections.
According to the actual parameters of the tunnel, the geometric model was established and the grid was divided. A numerical simulation of the wind field and pollutant concentration field in the tunnel in a turbulent state was conducted to establish a simulation model. The results showed that when the air flow moved to the front face of the tunnel, it produced a backflow to form an air flow vortex that rotated counterclockwise. When the carbon monoxide gas molecules were within the range of the turbulence, they moved irregularly with the wind field, which was not conducive to the diffusion of pollutants and was prone to the accumulation of polluting gases. The whole wind field tended to be stable at the plane 25 m from the tunnel's outlet. Pollutant gas diffused and migrated under the action of the wind field. Because the wind speed at the vent pipe orifice was the greatest and was the first to contact the pollutant gas, the concentration of gas was initially diluted. Turbulence in the wind field made it difficult for pollutants in the left-and right-hand corners of the roadway to diffuse.
According to this study, the results of the numerical CFD simulation model established by using the gas flow field data of the tunnel could be improved in the future to improve the ventilation facilities in the tunnel, and the measure of adding ventilation pipes are also to be preferred. According to the actual parameters of the tunnel, the geometric model was established and the grid was divided. A numerical simulation of the wind field and pollutant concentration field in the tunnel in a turbulent state was conducted to establish a simulation model. The results showed that when the air flow moved to the front face of the tunnel, it produced a backflow to form an air flow vortex that rotated counterclockwise. When the carbon monoxide gas molecules were within the range of the turbulence, they moved irregularly with the wind field, which was not conducive to the diffusion of pollutants and was prone to the accumulation of polluting gases. The whole wind field tended to be stable at the plane 25 m from the tunnel's outlet. Pollutant gas diffused and migrated under the action of the wind field. Because the wind speed at the vent pipe orifice was the greatest and was the first to contact the pollutant gas, the concentration of gas was initially diluted. Turbulence in the wind field made it difficult for pollutants in the left-and right-hand corners of the roadway to diffuse.
According to this study, the results of the numerical CFD simulation model established by using the gas flow field data of the tunnel could be improved in the future to improve the ventilation facilities in the tunnel, and the measure of adding ventilation pipes are also to be preferred.

Conclusions
This study simulated the concentration migration law of pollutants in the roadway. The main conclusions are as follows: (1) The simulated results of the wind field and pollutant concentration field in the actual roadway in a turbulence state produced by COMSOL software were as follows. When the air flow moved to the front face of the roadway, it was easy to generate a reverse flow to form an air flow vortex, where the air flow was stagnant. When the pollutant gas molecules came within range of the vortex, they rotated with the wind field, which was not conducive to the diffusion of pollutants, leading to the accumulation