A Mass-Conserving Model of Pollutant Transport in the Atmosphere Over Complex Terrain

A mathematical model and a verified numerical method are developed for the transport of passive and settling pollutants in the atmospheric boundary layer over complex terrain. The advection–diffusion equation is written in the terrain-following coordinate ???? = (????−ℎ)/(????????????????−ℎ) in fully conservative (Jacobian-weighted) form. The contravariant vertical velocity is diagnosed from the discrete continuity equation, so that the kinematic condition at the terrain surface is satisfied exactly and no mass can cross the ground. The horizontal wind, prescribed from a single-point forecast or reanalysis through a logarithmic profile in the height above ground, is made mass-consistent by a column-flux normalisation under a rigid lid. Turbulent diffusivities follow Monin–Obukhov similarity with a parabolic decay towards the top of the layer, and the stability parameter is obtained from the hourly solar radiation and air temperature. The finite-volume discretization uses shared face fluxes, first-order upwinding, and an explicit step chosen from a local positivity bound. It conserves mass to round-off, preserves a uniform field exactly, and keeps concentrations non-negative without clipping. Verification includes an analytical point-source solution (relative ????2 error 1.1%, ???? = 0.9999, ????????????2 = 1 on the finest grid), grid and time-step refinement, and a closed global budget with residual below 2×10−12. It is also shown that the widely used simplification of a horizontal wind (???? = 0) in the ????-equation violates the kinematic surface condition: over a 200 m ridge it removes more than 60% of the emitted mass within 15 min and produces a spurious lee-side attenuation of two to three orders of magnitude. With the conservative formulation, a ridge changes near-ground concentrations by 2–51%, the largest effect occurring for an elevated source under stable stratification, whereas atmospheric stability changes them by more than an order of magnitude. The model is coupled with the Copernicus digital elevation model and ERA5 reanalysis data and applied to the city of Tashkent.