#### Abstract

The underground ventilation in coal mine comes from the ground air so that the seasonal fluctuation of atmospheric temperature inevitably has influence upon its climate. The paper obtains a simplified mathematical calculation method of the radial and axial temperature fields to furnish a supplement to air temperature prediction of coal roadway with the emphasis on the influence of seasonal characteristic of inlet airflow. The results reveal that the propagation of radial and axial temperature wave has the similar characteristic of amplitude decay and phase delay. In the radial direction, the main factor is the dimensionless annual cycle time; and in the axial direction, it is influenced by the location-related dimensionless time, the Nusselt number, the Biot number, and the dimensionless annual cycle time. In addition, for the wave amplitude, they present the negative exponential trend as the location increases, the overall shape of which resembles a folded page; for the delayed phase, they present the linear variation trend as the location increases, the overall shape of which resembles a parallelogram. In the paper, the evolution mechanism of temperature field of coal roadway under seasonal fluctuation boundary is profoundly comprehended, which has important guiding significance for the field engineer.

#### 1. Introduction

Coal is the biggest energy sources in China. By statistics, it accounted for about 60.4 percent in the primary energy consumption in 2017 [1]. However, as the shallow coal resource is gradually reduced and even exhausted, the mining depth is growing year by year; and China's third resource prediction showed that the coal seam buried depth greater than 1000 m accounted for 59.5 percent of the total forecast. Therefore, it inevitably will suffer from more and more severe heat hazard menace in the not-too-distant future, and is directly linked to endanger the worker’s health, reduce the labor productivity, and increase the frequency of accidents and so on [2–4]. The situation has seriously restricted its safe and efficient exploitation so that government puts more efforts on the solution of the thermal damage problem.

The thermal damage governance of high-temperature mine mainly is based on two aspects [5]: the prediction of airflow temperature of coal mine and the implementation of reasonable cooling technical measures. Meanwhile, the accurate prediction for mine climate is a crucial factor that determines the optimal thermal control measures to evaluate the labor conditions and improve the mine environment. However, the variation rules of airflow temperature in underground mine are very complex because of the influence of many factors, principally including the heat dissipation from surrounding rock and the inlet initial temperature from atmosphere [6]. In this perspective, the research on the temperature characteristic under seasonal fluctuation boundary, which includes radial and axial temperature distributions of surrounding rock and airflow, has profound significance to perfect the basic theory of wind temperature prediction in high-temperature mine.

There have been many studies focused on the strata heat release and temperature field of underground coal mine roadway. With regard to the radial heat dissipation from surrounding rock, Щсрбаиъ. АП firstly [7] proposed a dimensionless parameter to quantify the heat dissipating capacity. Then an analytical solution for the parameter was established by Cen et al. [8] and he pointed that it was composed of Bessel function and Kelvin function. Sun [9] proposed a simplified analytic expression to calculate the parameter. Other studies also identified that it was an important physical quantity which described the heat dissipation of roadway surrounding rock [10, 11]. After that, Qin et al. [12–14] presented many relation curves between the parameter and Fourier modulus and analyzed the strata heat release using numerical simulation methods. Zhang et al. [15] established a theoretical model of thermal diffusion using finite different method, obtaining the temperature distribution of surrounding rock under the coupling effect of airflow and seepage. At the same time, in the application of business software, Bai et al. [16] analyzed the heat transfer between surrounding rock and airflow based on Fluent software in the case of dry roadway. Danko et al. [17] investigated the transfer phenomena of heat and mass on the wall of roadway based on CLIMSIM software. In the experimental aspect, Zhang et al. [18–20] deduced that the similarity criterion of heat conduction of roadway surrounding rock was the Fo (Fourier number), the Bi (Biot number), and the R (Radius of roadway) in the heat transport process and the Nu (Nusselt number), the Re (Reynold number), and the Pr (Prandtl number) in the heat convection process between surrounding rock and airflow and then certified the feasibility of modeling experiment method. Moreover, the influence of specific factors to the characteristic of temperature field also was analyzed, including the virgin rock temperature [21], the moisture of roadway wall [22–25], the section shape of roadway [26], heterogeneous rocks and anisotropy [27], and ventilation time and ventilation volume [28]. These above factors reflected their influences to heat release of surrounding rock from different sides, which involved the disturbance range of airflow temperature, heat dissipating capacity of surrounding rock, and the distribution rule of temperature field.

In addition, for the axial temperature field, Roy et al. [29] developed a computer program to predict the time-dependent climatic change in an underground airway. After that, Liu [30] presented a method of approximate calculation of airflow temperature in mine based on the heat radiation from country rock, the depth of mine and the humidity of the airflow. Gao et al. [31] analyzed the influence of water evaporation place distribution on temperature and humidity of airflow. Xie [32] discussed some influence factors through the accurate observation and analysis of the thermodynamic parameters. Liu [33] deduced the basic law of temperature change in an intake shaft and pointed that it changed with the temperature of airflow at the inlet and the shaft depth. Roy [34] presented a simplified computational approach for prediction of transient climatic conditions in underground mine airways.

It was well known that there were two aspects of the periodic change of atmospheric temperature: diurnal variation in a day and annual variation in a year. The ventilation of coal mine came from the ground air so that the regular fluctuation of atmospheric temperature predicatively had influence on its climate. The variation time of diurnal periodicity was relatively short compared with the annual periodicity, and its temperature difference was also small, so it had little influence on the airflow temperature in mine; by contrast, the field measurement indicates that the mining face had obvious annual periodic fluctuation characteristics [35]. However, as mentioned above, most of research only took into account the properties of surrounding rock media and took the temperature of airflow as a constant but ignored its seasonal fluctuation factor. Yi et al. [36] proposed that cooling measures of high-temperature mine should take different implementation on the basis of the heat hazard severity of the season, but without revealing when to implement the optimal cooling methods.

In this paper, a technique of variable separation method is introduced to decompose the air inlet temperature of coal mine into two portions: fluctuation boundary and steady boundary. Here, as to the analysis above, we only discuss the influence of seasonal fluctuation of inlet airflow temperature to the radial and axial temperature fields of coal roadway. Through comparing with the analytic solution of a half-infinite plane on even ground, we previously [37] analyzed the mathematical solution of radial temperature field of surrounding rock in detail, and further take a supplement to it, the correctness of which is verified by the experiment in the paper. On this basis, the mathematical solution of axial temperature distribution can be obtained from the law of conservation of energy. Ultimately, all the factors that affect the wave amplitude and phase angle of inlet airflow temperature are synthetically evaluated from the whole space-time perspective when to implement the suitable cooling measures.

#### 2. Mathematical Model

##### 2.1. Fundamental Assumptions

In the present work, the following basic assumptions have been made during the implementation of heat dissipation model of surrounding rock.

(i) As shown in Figure 1, the cross-section of roadway is simplified as the shape of circular and the heat transmission within the surrounding rock is only along radial direction. Γ_{1} and Γ_{2} represent the inner boundary and the outer boundary, respectively.

(ii) The surrounding rockmass is homogeneous and isotropic, the original temperature field of which is a constant.

(iii) The effect of moisture content and thermal radiation of surrounding rock wall is ignored in this study.

(iv) The airflow is laminar and incompressible, and its pressure does not vary as temperature; the air inlet temperature of coal mine is considered as a trigonometric function, which is shown in Figure 2.

##### 2.2. Governing Equation

The governing equation of transient heat conduction of roadway surrounding rock can be established based on the law of conservation of energy and the Fourier's law in the cylindrical coordinate.

The boundary conditions are

The initial conditions are

where is the temperature field of surrounding rockmass, °C; is the thermal conductivity of rock, W/(m•°C); is the density of rock, Kg/m^{3}; is the specific heat capacity of rock, J/(Kg*∙*°C); is the initial temperature of surrounding rockmass, °C; r is the radial depth, m; h is the convective heat transfer coefficient, W/(m^{2}•°C); q is the heat flux, W/m^{2}; is the seasonal airflow temperature, °C; is the average annual temperature, °C; t is the time of heat exchange between airflow and surrounding rockmass, s; *τ* is the time of one year (365×24×3600 s). is the amplitude of airflow temperature; *φ* is the phase.

##### 2.3. Dimensionless Analysis

For simplicity, the following dimensionless numbers can be introduced to make it more widely applied:

where *θ* is the dimensionless excess temperature; R is the dimensionless radius; r_{0} is the radius of roadway, m; is the thermal diffusivity of surrounding rock, m^{2}/s; Bi is the Biot number; Fo is the Fourier number; is the dimensionless annual cycle time; is the dimensionless fluctuation amplitude.

Therefore, (1) can be simplified as

The boundary and initial condition also can be similarly obtained based on the above analysis, which can be expressed as

##### 2.4. Discretized the Governing Equation

Obviously, it is very important that the property of boundary condition is properly analyzed to obtain the accurate result. It is well known that there are three common boundary conditions for heat conduction problem, the first-, second-, and third-type boundary conditions (i.e., Dirichlet, Neumann and Robin, respectively). Compared with the last one, the first two is easily intelligible, which can be taken as its particular case. Unfortunately, the seasonal fluctuation boundary belongs to a more special case. The heat dissipation problem under this condition becomes so complex that it is difficult to implement the numerical simulation straightway. On this ground, we introduce a technique of variable separation method to decompose the complicated fluctuant boundary into two simply parts: one for the average annual temperature of airflow and the other for the periodic temperature with characteristic of cosine wave . In addition, this cosine wave can also be divided into two aspects: one for the wave amplitude and the other for the cosine function . Note that the physical quantity is a nondimension parameter, which can be defined as

Similarly, (1) also consists of two parts: one represents the influence of the average annual temperature to the temperature filed of surrounding rock; the other represents the influence of the periodic temperature . By that analogy, (5) can also be divided into two parts: and .

where

Therefore, the corresponding definite condition of and can be defined as

So the problem for (1) is converted to solve the formulation of (10) and (11), respectively. To the best knowledge of the authors, most of the research have studied the impact of the average annual temperature but neglecting the periodic fluctuation factors . For this reason, we only discuss the variation regularity of temperature field of roadway surrounding rock under the influence of seasonal fluctuation factors in this article.

#### 3. Numerical Analysis

##### 3.1. Physical Model Discretization

The given physical model of temperature filed is simplified as a circular roadway, as shown in Figure 1, the problem of which has been transformed to solve the axisymmetric heat conduction. Therefore, the discretization of mesh for the physical model can be shown in Figure 3, and the calculated domain from 0 to N depicts from the wall of roadway to a certain depth of surrounding rock. The bold line near the centre denotes the coal roadway, and each thin line represents every cutoff point of the gridding and is divided into N sections (i.e., N=1,2,3,4…). The rule of the gap distance between two adjacent gridding nodes obeys the geometric progression due to the intense heat exchange near the roadway boundary. At the same time, the imaginary line in Figure 3 is located in the middle of two adjacent gridding nodes, and any two adjacent imaginary lines represent one integral domain, which can be considered as a complete cylinder ring. But it can be seen that the integral domain model of the boundary for the node 0 and the node N, comparing with the internal, should be defined as half of the cylinder ring since it only contains one imaginary line.

##### 3.2. Mathematical Model Discretization

The thermal balance equations of each cylinder ring within the surrounding rock can be established based on the law of energy conservation, which can be expressed as follows:

In the same way, the node 0 and the node N can be displayed as well:

where the subscripts i-1, i, and i+1 denote the node number of the corresponding position; the superscripts k and k-1 denote the time nodes; ΔFo denotes the time step at the .

Therefore, a complete set of heat conduction equation of surrounding rock has been established by three equations, (12), (13), and (14).

The numerical calculation program is developed using VB software based on the aforementioned analysis. In order to reduce the error, the independence of grid and time step has been verified. At last, the radial grid based on the equal proportion is divided by 50 nodes, and the proportionality coefficient is 1.1. The time step is divided as the equal interval, the spacing of which is 1000^{−1} of the dimensionless annual cycle.

#### 4. Result

##### 4.1. Mathematical Model of Temperature Distribution along the Radial Direction

To study the regularity of periodic temperature fluctuation of roadway surrounding rock mass, we refer to the temperature distribution function of a semi-infinite medium on even ground [38].

where

From the analysis of (15), the amplitude *α* and phase* η* is closely related to Bi and Fo. Moreover, the attenuation rate of amplitude

*α*follows the law of negative exponent as the position X increases and the retardation rate of phase

*is linearly related to the position X. The correlation coefficients between them are the function of wave period . It can be seen that the temperature field of surrounding rock under periodic fluctuation has the similar characteristics to the semi-infinite medium’s. Therefore, assume the nondimensional temperature field of surrounding rock under periodic boundary has the similar mathematical expressions, which can be defined as*

*η*where X represents the depth from the surface of surrounding rock to the interior and X=R-1.

Based on (17), the wall temperature wave (X=0) can be obtained as

By analyzing and comparing with the numerical result of (11), we obtain the similar functional relations between A, B, C, D, and , Bi for (17), which can be expressed using (19). The errors of them relative to the theoretical results are all less than 1%. The detailed derivation process can be referred to in the literature [37].

##### 4.2. Experimental Verification

The field test in coal mine has certain difficulties to verify the correctness of numerical calculation results due to the limitation of objective conditions like expensive cost, long period, hostile site environment, and so on. Therefore, a method of similarity experiment has been implemented. In addition, because the axial temperature field is the result of theoretical derivation from the radial temperature distribution based on the above analysis, here, we will only verify the temperature distribution rule of the radial direction.

###### 4.2.1. Experimental System

As shown in Figures 4 and 5, the experimental platform which is designed according to the similarity criterion is composed by three systems, experimental mimic system, air temperature control system, and temperature field monitoring system, respectively. The geometric parameter of mimic system welded with stainless steel plates is 1000 mm ×1000 mm ×500 mm, and the radius of analogous roadway is 50 mm; the outer surface of test model is covered by the thermal insulation material. The design parameter of control system for inlet air temperature involves the control range -20~150°C, deviation ≤ ±1°C, accuracy ± 0.1°C and air velocity 0~6 m/s. The temperature monitoring system includes the temperature sensor and the data acquisition unit Datataker800, which can be implemented by collecting temperature signal in real time. Figure 6 shows the measuring points layout of thermal resistance in the experimental mimic system, including one in the center of roadway to control and monitor the inlet air temperature, one at the surface of roadway to measure the wall temperature, and the rest within the surrounding rock, the overall outline of which likes a spiral in order to reduce the interplay between sensors.

In the process of experiment, the experimental installation cannot keep the inlet air temperature change continuity due to the limitation of its own performance. So in order to reduce its work loads and make the inlet air temperature accord with the nature of the trigonometrical functions, the implement of temperature control in each cycle is averagely divided into 24 pieces based on (3), which then is inputted into the temperature control system.

In addition, since the experiment period is relatively short and the temperature of indoor is less affected by that of the outside, we assume that the temperature within surrounding rock model is equal to that of the ambient in initial moment, which can be approximated as a constant 25°C. Thus, the average temperature of inlet air is also set as 25°C and its impact can be ignored. In this paper, the plaster is selected as the main material of similar simulation experiment, in which the powder of iron oxide is added to reinforce its heat conductivity, which is shown in Figure 7. The wave amplitude is set as 20°C.

**(a) Gypsum**

**(b) Iron oxide**

Then, Table 1 lists the thermal physical parameters of the test model, and the experimental schemes can be shown in Table 2. Moreover, the convective heat transfer coefficient between air and surrounding rock model can be calculated through (20), which is based on the Dittus-Boelter formula [39].

where V_{f} is the velocity of airflow, m/s; is the Kinematic viscosity of airflow, m^{2}/s.

###### 4.2.2. Mathematical Model Validation

The wave amplitude and delay phase of numerical calculation are compared with the experimental model. It can be observed from Figures 8–11 that the computational results of mathematical model match with the experimental results very well, which shows that this calculation method has high computational accuracy in this paper.

It can be seen from Figures 8 and 9 that the amplitudes of temperature wave show the law of negative exponent as the radial depth increases, and from Figures 10 and 11 that the phases show the linear variation.

Furthermore, Figures 8 and 10 present that the amplitude and phase of the wall temperature wave (R=0) are very little influenced by the cycle; but in the interior (R>0), the shorter the period, the greater the rate of decay and delay. It is well-know that the mechanism of heat transfers that the high-temperature object transmits its internal energy to the low-temperature object. Hence, when the period of temperature fluctuation has been shortened, the temperature difference per unit time has been increased so that the rate of heat transfer also increases, the characteristic of which exists in the wall and interior. So it not only affects the wall temperature wave of surrounding rock but also the internal.

In addition, Figures 9 and 11 present, for the wall temperature wave(R=0), that the greater the wind speed, the larger the amplitude but the smaller the phase. Meanwhile, within the surrounding rock (R>0), the gradient of the decay or delay under the condition of different wind speeds is exactly the same. Therefore, the wind speed only affects the wall temperature of surrounding rock but not the internal. The reason is that the heat transport mode between the wall of surrounding rock and air is through convection, which has greatly related to the speed of air current from (20); that is, the parameter A has a close relationship with Fo, Bi; but in the interior, it is through the conduction, which has only relationship with the thermal diffusivity of surrounding rock, as well the parameter C is only a function of Fo.

##### 4.3. Mathematical Model of Temperature Distribution along the Axial Direction

According to the above-mentioned analysis, the interior of surrounding rockmass is going to form a temperature field that fluctuates periodically after the seasonal airflow ventilating into the underground mine, and in turn the heat release from the surrounding rock is bound to affect periodically the air. Therefore, the airflow temperature distribution along the axial direction is the result of the interaction between the inlet airflow and the heat release from rock, which makes the air temperature prediction complicated.

Based on the theoretical analysis and field observation, we can assume the axial temperature distribution as follows:

where P is the amplitude; Q is the phase. It is obvious that they both are associated with the concrete position in the roadway and not related to the time when the inlet air temperature is determined. At the same time, the wall temperature of the same location can be obtained from (18), which can be expressed as follows:

In the absence of moisture evaporation, all of the heat lost from the surrounding rockmass is gained by the wind. Therefore, as shown in Figure 12, a heat-balance equation in the unit length of dL can be established.

where G is the mass flow of air, m^{3}/s; C_{f} is the specific heat of air, J/(Kg•°C); U is the section perimeter of roadway, m; dT_{af} is the temperature difference in the distance of dL; L is the length from the ventilation inlet to a site used the airflow, m.

where S is the section area of roadway, m^{2}; *ρ*_{f} is the density of air, Kg/m^{3}.

Substituting (21), (22), and (2) into (23), we can obtain that

where

Simplify the above equation

It should be noticed that (27) is strict equality in both side at any time and any location. At the same time, there is no correlation between the amplitude P and the phase Q. Therefore, the following formula can be obtained:

Further simplification, we can obtain that

where P_{0} and Q_{0} are the initial amplitude and the phase in the inlet of roadway, respectively. From (3), that is to say,

From the discussion above, the calculation formula of amplitude and phase of the air temperature along the axial direction has been deduced. So the conclusion can be obtained that the amplitude shows the law of negative exponent along the axial direction of the roadway, and the phase changes linearly.

Substituting (29) and (30) into (21), the axial temperature distribution of roadway can be obtained as follows:

because of

where represents the Nusselt number; represents the dimensionless cycle time that related location. Therefore, in the axial direction, there also is a dimensionless cycle time on the synchronous occurrence of the maximum value of temperature, and it can be calculated as follows based on (31) and (32):

As a result, we can obtain the dimensionless expression of temperature distribution along the roadway.

##### 4.4. Project Application of Mathematical Model of Temperature Field

###### 4.4.1. Calculation Scheme

This paper takes a specific mine as an example in order to systematically analyze the influencing factors of temperature field. The relevant parameters are shown in Table 3 [40].

However, for any specific coal mine, the main factors that can be controlled by human are generally only the ventilation velocity and the section area of roadway, the influence of which on the amplitude and phase is analyzed based on the coal mine safety regulations. It stipulates that the wind speed in the main intake and return airways should not exceed 8 m/s, and the net section of roadway must meet the needs of pedestrian, transportation, ventilation, safety facilities, equipment installation, overhaul, and construction. Therefore, the following calculation schemes can be listed in Table 4.

###### 4.4.2. Result Analysis

It can be seen from Figures 13 and 14 that the farther the entrance of roadway along the axial direction, the larger the amplitude of temperature wave, and it shows the law of negative exponent as the axial depth increases. In addition, the cumulative distance from the ventilation inlet to the working face is generally approximately 6000 m, and the variation range of the amplitude is around 2~6°C, which is consistent with that of field observation at the working face [41].

It can be seen from Figures 15 and 16 that the farther the entrance of roadway along the axial direction, the more delayed the phase of temperature wave, and it shows the linear variation rule as the axial depth increases. In addition, when the cumulative distance is about 6000m, the variation range of the delayed phase is around 0.5~2, namely, the temperature wave of the working face lag behind that of the ventilation inlet about 30~116 days. At the same time, the field observations at the working face shows that the maximum value of ground temperature in a year appears in the mid-to-late July and that of working face in early September [38], which agrees with the calculation.

Figures 13 and 15 present that the bigger the wind speed, the smaller the attenuation rate of wave amplitude and the slower the delay rate of the phase. The reason is that the heat exchange time between air and surrounding rock is shortened when the wind speed increased, so the attenuation rate and delay rate at unit time is decreased.

Figures 14 and 16 present that the bigger the radius of roadway, the smaller the attenuation of wave amplitude and the slower the delay rate of the phase. Due to the fact that the wind speed is a constant, the rate of heat exchange is also assumed to be constant. Moreover, the ventilation quantity in unit time inevitably increases as the roadway radius enlarges, which lead to a decreased attenuation rate and delay rate at unit length.

On the incorporation of the radial and axial temperature distribution, the synthesis analysis of temperature field in three-dimensional space can be obtained when the wind speed is 3m/s and the radius is 3m.

Figure 17 shows that the surrounding rock and air have the attenuation effect on the amplitude of temperature wave, and both present the law of negative exponent as the corresponding location increases, the overall shape of which resembles a folded page. Figure 18 shows that they both have the postponed effect on the phase and present the linear variation trend as the corresponding location increases, the overall shape of which resembles a parallelogram.

Moreover, compared with the air, the decay rate or delay rate of temperature wave within surrounding rock is relatively fast, which is affected by the dimensionless annual cycle time (thermal diffusivity of rock, radius of roadway, cycle time, and radial depth). However, in the axial direction, they are influenced by the location-related dimensionless time (thermal diffusivity of air, radius of roadway, and air velocity), the Nusselt number, the Biot number, and the dimensionless annual cycle time. Therefore, whether the decay of amplitude or the delay of phase, the effect of the rock is greater than that of the airflow.

Figure 19 shows that the ventilation distance has great influence on the amplitude and the phase of air temperature wave. The location-related cycle time along the axial direction can be calculated based on (33) and has close connection with the ventilation distance and velocity. In this instance, it is about 4.3 hours or the distance is 46600m. As mentioned previously, the ventilation rate in mine is relatively small and the maximum is 8m/s. Therefore, it is difficult to achieve a complete cycle fluctuation along the axial direction at one time.

From Figure 20, the temperature of air current may be greater than that of surrounding rock in the hottest season, and then the air could dissipate constantly heat to the surrounding rock. The other way round, the coldest season could present the opposite characteristic. That how to make good use of these features to govern the mine heat damage needs further discussion.

#### 5. Conclusions

In this work, the evolution mechanism of temperature field of coal roadway under seasonal fluctuation boundary is elaborated in detail, involving the radial and axial direction, to better implement the air temperature prediction in the underground high-temperature mine. The main conclusions can be drawn.

The seasonal characteristics of the inlet air temperature can be regarded as a simple harmonic wave, and one dimensional unsteady heat conduction differential equation and its periodic boundary condition is established. The heat-balance equation is derived by finite volume method.

The correctness of the temperature distribution function of surrounding rock of roadway is verified by the laboratory experiments. In addition, for the wall temperature wave, the greater wind speed, the larger wave amplitude but the smaller delayed phase; and they are almost not influenced by the cycle period. But in the interior, they appear to be the same attenuation or delay trend under the conditions of different wind speed; and at the same time the shorter the period, the greater the rate of decay and delay.

On this basis, the mathematical model of temperature distribution along the axial direction is obtained from the heat-balance equation. The wind speed and roadway radius have great influence on the amplitude and the phase of temperature wave of air. The bigger the wind speed, the smaller the attenuation rate of wave amplitude and the slower the delay rate of phase. The bigger the roadway radius, the smaller the attenuation rate of wave amplitude and the slower the delay rate of phase.

The radial and axial propagations of temperature wave have the similar characteristics of the amplitude decay and phase delay. For the wave amplitude, they both present the law of negative exponent as the location increases, the overall shape of which resembles a folded page. For the phase angle, they both present the linear variation trend as the location increases, the overall shape of which resembles a parallelogram.

Compared with the air, the decay rate or delay rate of temperature wave within surrounding rock is relatively fast, which is affected by the dimensionless annual cycle time (thermal diffusivity of rock, radius of roadway, cycle time, and radial depth). However, in the axial direction, they are influenced by the location-related dimensionless time (thermal diffusivity of air, radius of roadway and air velocity), the Nusselt number, the Biot number, and the dimensionless annual cycle time. Therefore, whether the decay of amplitude or the delay of phase, the effect of the rock is greater than that of the airflow.

In the paper, the seasonal variation of the airflow temperature in mine causes the problem of seasonal heat disaster; the principle analyzed scientifically can be guided engineering practice.

#### Data Availability

The data used to support the findings of this study are included within the supplementary information file.

#### Conflicts of Interest

There are no conflicts of interests for all authors of this paper.

#### Acknowledgments

This paper was supported by the University Key Scientific Research Programs of Henan Province (19B620002), the Key Research and Development Plan’s Special Projects of Henan province (182102310772), Henan Basic and Frontier Technology Research Project (162300410070), and the Doctoral Program of Zhengzhou University of Light Industry (2016BSJJ048).

#### Supplementary Materials

The supplementary material consists of three files, one for related data and other two for related graphics.* (Supplementary Materials)*