Abstract: Heat dissipation from the surrounding rock is one of the major sources of thermal hazards in excavation roadways. Accurate determination of the axisymmetric temperature field of the surrounding rock is therefore essential for quantifying rock heat dissipation and predicting airflow temperature. To optimize the control volume division scheme in the finite volume method for this axisymmetric heat dissipation problem, software was developed on the Visual Studio platform to compute the dimensionless temperature field in excavation roadways. Based on this software, the evolution of the temperature field, the unsteady heat transfer number under different division schemes, and the influence of grid shape on the numerical results were systematically investigated, and an appropriate division scheme was identified. The results show that the cooling zone is relatively limited near the heading face, but becomes much larger in regions farther away. Under the same conditions and at the same location, the cooling radius predicted by the 1/3 division scheme is the largest, followed by the 1/2 division scheme, while the 2/3 and center of gravity division schemes give slightly smaller values with little difference between them. A similar trend is observed for the temperature field distribution. The unsteady heat transfer number reaches its maximum near the heading face in the peripheral roadway and decreases gradually with increasing distance from the heading face. The values obtained using the 2/3 and center of gravity division schemes are the largest and remain close to each other, whereas the 1/2 division scheme yields the smallest values. For the center of gravity division scheme, the difference in the unsteady heat transfer number between the coarse and refined grids at the heading face is only 0.1269, which is the smallest among all schemes, indicating better grid adaptability. In addition, the center of gravity division scheme fully covers the geometric domain without omission, overlap, or extension beyond the computational region, further confirming its rationality. The simulated results agree well with the measurements in terms of variation trend, showing that the surrounding rock temperature increases with borehole depth and eventually approaches a stable value. For the two boreholes, the maximum relative errors between simulated and measured temperatures are 12.71% and 13.09%, respectively, while the corresponding average relative errors are 3.84% and 3.94%. These results satisfy engineering requirements and demonstrate the validity of both the thermal field mathematical model and the discrete model based on the center of gravity division scheme.