## 2 Mathematical Model  
### Governing Equations and Numerical Algorithm  
Previous studies have shown that the standard k-\(\varepsilon\) turbulence model remains the most widely used approach to investigate the dispersion characteristics of vehicular pollutants and model atmospheric dispersion problems [18, 19, 20, 21]. Therefore, the Reynolds-averaged Navier-Stokes equation (RANS) [22] and standard k-\(\varepsilon\) models were adopted, and the following basic principles and assumptions were followed in the calculation:  
1. Because the pollutant discharge rate was far lower than the speed of sound at the upper vents, the compressibility of discharge airflow can be ignored. This process can be regarded as the steady turbulent flow of incompressible fluid [23]. In addition, the influence of ambient temperature on the discharge airflow was temporarily omitted. 2. The fluid was set as a continuum, and no spaces existed between particles. The total fluid mass flowing through each cross section per unit time remained unchanged. 3. Considering the higher contribution rate of NO\({}_{2}\) than of CO, NO\({}_{2}\) was used as the main source of pollution to study the diffusion of pollutants. Based on the rationale and assumption reported above, the equations for controlling the pollutant diffusion at the upper vent of urban tunnels considering the buoyancy effect of turbulent jets are as follows:  
Continuity equation:  
\[\frac{\partial}{\partial x_{i}}(\rho u_{i})=0, \tag{1}\]  
Momentum conservation equation:  
\[\frac{\partial}{\partial x_{j}}(\rho u_{i}u_{j})=-\frac{\partial p}{\partial x_ {i}}+\frac{\partial}{\partial x_{j}}\left[\mu\frac{\partial u_{i}}{\partial x _{j}}-\rho\overline{u_{i}^{\prime}u_{j}^{\prime}}\right]+S_{i} \tag{2}\]  
\(k\) equation:  
\[\frac{\partial}{\partial x_{i}}(\rho ku_{i})=\frac{\partial}{\partial x_{j}} \left[\left(\mu+\frac{\mu_{t}}{\sigma_{k}}\right)\frac{\partial k}{\partial x _{j}}\right]+G_{k}+G_{\text{b}}-\rho\epsilon \tag{3}\]  
\(\epsilon\) equation:  
\[\frac{\partial}{\partial x_{i}}(\rho\epsilon u_{i})=\frac{\partial}{\partial x _{j}}\left[\left(\mu+\frac{\mu_{t}}{\sigma_{\epsilon}}\right)\frac{\partial \epsilon}{\partial x_{j}}\right]+C_{1\epsilon}\frac{\epsilon}{k}(G_{k}+C_{3 \epsilon}G_{b})-C_{2\epsilon}\rho\frac{\epsilon^{2}}{k} \tag{4}\]  
Species mass conservation equation:  
\[\frac{\partial(\rho c_{s}u_{i})}{\partial x_{i}}=\frac{\partial}{\partial x_ {j}}\left[D_{s}\frac{\partial(\rho c_{s})}{\partial x_{j}}\right]+S_{s} \tag{5}\]  
where \(G_{k}=\mu_{t}(\partial u_{i}/\partial x_{j}+\partial u_{j}/\partial x_{i}) \partial u_{i}/\partial x_{j}\) is the generation term of turbulent kinetic energy;\(\mu_{t}=\rho C_{\mu}k^{2}/\epsilon\) is the turbulent viscosity coefficient; \(G_{b}=\left(-\mu_{t}g_{i}/\sigma_{\rho}\right)\partial\rho/\partial x_{i}\) refers to the buoyancy correction term. The variables in the above formula are defined in Table 1. The model coefficients \(C_{\mu}\), \(\sigma_{k}\), \(\sigma_{\epsilon}\), \(C_{1\epsilon}\), \(C_{2\epsilon}\), and \(C_{3\epsilon}\) were determined according to the principle of gradualism, and their values are shown in Table 2 [24]. The above governing equations were discretized in Ansys Fluent 2019R3 by using the finite volume method in a fully implicit way. The convective term adopted the second-order