Chapter 16.2: Advection-Diffusion & OSPM

16.2.1 The Advection-Diffusion Equation

The concentration $C(\mathbf{x},t)$ of a passive scalar (e.g. NO₂, PM₂.₅) in a street canyon is governed by the advection-diffusion equation with a source term:

$$\frac{\partial C}{\partial t} + \mathbf{u} \cdot \nabla C = \nabla \cdot (K \nabla C) + S$$

Each term has a physical interpretation:

$\partial C / \partial t$: local rate of change of concentration
$\mathbf{u} \cdot \nabla C$: advection by the mean wind field
$\nabla \cdot (K \nabla C)$: turbulent diffusion with eddy diffusivity $K$
$S$: source term from vehicle emissions at street level

This is the Reynolds-averaged form, where $K$ encapsulates the effect of turbulent fluctuations on scalar transport. In a street canyon, $K$ varies strongly with position, especially with height.

16.2.2 Eddy Diffusivity from Mixing Length Theory

Prandtl’s mixing length theory gives the turbulent viscosity as $\nu_t = \ell_m^2 |dU/dz|$, where the mixing length near a wall is $\ell_m = \kappa z$. Inside the canyon, the flow is bounded both by the ground ($z = 0$) and the roof level ($z = H$). The mixing length must vanish at both boundaries, giving:

$$\ell_m(z) = \kappa \sqrt{z(H - z)}$$

The eddy diffusivity for a scalar is related to the turbulent viscosity through the turbulent Schmidt number $Sc_t \approx 0.7$:

$$K(z) = \frac{\nu_t}{Sc_t} = \frac{\kappa\,u_*\,z\,(1 - z/H)}{Sc_t}$$

where $u_*$ is the friction velocity associated with the canyon vortex. This parabolic profile has key features:

Maximum diffusivity at mid-height: $K_{\max} = \kappa u_* H / (4\,Sc_t)$
Zero at ground and rooftop—explaining why pollutants are trapped near both surfaces
The friction velocity $u_*$ connects to Module 16.1 through $u_* \sim u_v$

16.2.3 Gaussian Plume Green’s Function

For steady-state conditions with constant wind speed $u$, the Green’s function of the advection-diffusion equation in unbounded space is the Gaussian plume:

$$C(x,y,z) = \frac{Q}{2\pi\,\sigma_y\,\sigma_z\,u} \exp\!\left(-\frac{y^2}{2\sigma_y^2}\right) \exp\!\left(-\frac{(z - z_s)^2}{2\sigma_z^2}\right)$$

where $Q$ is the emission rate (g/s), $z_s$ is the source height, and the dispersion parameters grow with downwind distance:

$$\sigma_y(x) = a_y\,x^{b_y}, \qquad \sigma_z(x) = a_z\,x^{b_z}$$

Method of Images for Bounded Canyons

In a street canyon, the ground ($z = 0$) and the roof level ($z = H$) act as reflecting boundaries. The no-flux condition $\partial C / \partial z = 0$ at both surfaces is satisfied by the method of images: we place an infinite series of image sources at $z = \pm 2nH \pm z_s$:

$$C_{\text{canyon}} = \frac{Q}{2\pi\sigma_y\sigma_z u} e^{-y^2/(2\sigma_y^2)} \sum_{n=-\infty}^{\infty} \left[ e^{-(z - z_s - 2nH)^2/(2\sigma_z^2)} + e^{-(z + z_s - 2nH)^2/(2\sigma_z^2)} \right]$$

In practice, the series converges rapidly; terms with $|n| \leq 3$ are sufficient for $\sigma_z < H$.

16.2.4 The Operational Street Pollution Model (OSPM)

The OSPM, developed by the Danish National Environmental Research Institute, decomposes the street-level concentration into three additive contributions:

$$C_{\text{street}} = C_D + C_R + C_B$$

Direct Contribution C_D

The direct contribution is from the Gaussian plume of emissions carried by the along-canyon wind component. For a line source of strength $q_L$ (g/m/s):

$$C_D = \frac{q_L}{\sqrt{2\pi}\,\sigma_w\,u_s} \cdot f_d$$

where $\sigma_w$ is the vertical velocity fluctuation, $u_s$ is the street-level wind speed, and $f_d$ is a geometric factor.

Recirculation Contribution C_R

The recirculation zone traps pollutants behind the upwind building. OSPM models this as:

$$C_R = \frac{q_L \cdot L_R}{u_v \cdot W \cdot H} \cdot f_r$$

where $L_R \approx \min(W, H)$ is the recirculation length and $f_r$ accounts for the residence time.

The Trapping Ratio

The ratio of recirculation to direct contribution reveals why deep canyons are dangerous:

$$\frac{C_R}{C_D} \sim \frac{H}{W} = \alpha$$

As the canyon deepens ($\alpha$ increases), the recirculation contribution grows linearly. For a typical Parisian boulevard with $\alpha \approx 1$, recirculation doubles the street-level concentration. For narrow streets in old European cities with $\alpha \approx 3$, the recirculation contribution dominates by a factor of three.

Background Concentration C_B

The background term $C_B$ represents the urban background level measured at rooftop stations, typically 20–40 μg/m³ for NO₂ in European cities. This is added as a constant.

16.2.5 Python: OSPM Concentration vs Canyon Geometry

This code implements the OSPM decomposition and plots how each contribution changes with canyon aspect ratio.

OSPM Concentration vs Canyon Geometry

Python

script.py118 lines

import numpy as np
import matplotlib.pyplot as plt

# ─── OSPM Parameters ───
q_L = 0.001          # line source strength (g/m/s) — typical for moderate traffic
sigma_w = 0.8        # vertical velocity fluctuation (m/s)
u_H = 5.0            # rooftop wind speed (m/s)
kappa = 0.41         # von Karman constant
C_D_drag = 0.005     # drag coefficient for vortex velocity
C_B = 25.0           # background NO2 (ug/m3)
W = 20.0             # canyon width (m)

# ─── Aspect ratio sweep ───
alpha_arr = np.linspace(0.1, 4.0, 200)
H_arr = alpha_arr * W

# Vortex velocity (from Chapter 16.1)
u_v_arr = C_D_drag * u_H / ((1 + alpha_arr) * kappa)

# Street-level wind speed (fraction of rooftop)
u_s_arr = np.maximum(u_v_arr * 0.5, 0.3)  # min 0.3 m/s for calm conditions

# ─── OSPM Components ───
# Direct contribution (Gaussian plume)
f_d = 1.0  # geometric factor ~ 1
C_direct = (q_L / (np.sqrt(2 * np.pi) * sigma_w * u_s_arr)) * f_d * 1e6  # convert to ug/m3

# Recirculation contribution
L_R = np.minimum(W, H_arr)  # recirculation length
f_r = 2.0  # residence time factor
C_recirc = (q_L * L_R / (u_v_arr * W * np.maximum(H_arr, 1.0))) * f_r * 1e6

# Total
C_total = C_direct + C_recirc + C_B

# Trapping ratio
trap_ratio = C_recirc / C_direct

# ─── Plot 1: OSPM Decomposition ───
fig, axes = plt.subplots(1, 3, figsize=(18, 5))

ax1 = axes[0]
ax1.fill_between(alpha_arr, 0, C_B * np.ones_like(alpha_arr),
                 color="#064e3b", alpha=0.5, label="$C_B$ (background)")
ax1.fill_between(alpha_arr, C_B, C_B + C_direct,
                 color="#10b981", alpha=0.6, label="$C_D$ (direct)")
ax1.fill_between(alpha_arr, C_B + C_direct, C_total,
                 color="#f59e0b", alpha=0.6, label="$C_R$ (recirculation)")
ax1.plot(alpha_arr, C_total, color="white", linewidth=2, label="$C_{total}$")

# WHO limit
ax1.axhline(40, color="#ef4444", linestyle="--", linewidth=1.5, label="WHO NO$_2$ limit (40 μg/m³)")

ax1.set_xlabel("Aspect Ratio α = H/W", color="white", fontsize=12)
ax1.set_ylabel("NO₂ Concentration (μg/m³)", color="white", fontsize=12)
ax1.set_title("OSPM Concentration Decomposition", color="white", fontsize=13, fontweight="bold")
ax1.legend(fontsize=9, facecolor="#1e293b", edgecolor="#334155", labelcolor="white", loc="upper left")
ax1.set_facecolor("#0f172a")
ax1.tick_params(colors="gray")
ax1.grid(True, alpha=0.2, color="#334155")
for spine in ax1.spines.values():
    spine.set_color("#334155")

# ─── Plot 2: Trapping Ratio ───
ax2 = axes[1]
ax2.plot(alpha_arr, trap_ratio, color="#06b6d4", linewidth=2.5)
ax2.plot(alpha_arr, alpha_arr, color="#94a3b8", linewidth=1, linestyle="--", label="$C_R/C_D = \\alpha$ (theory)")
ax2.set_xlabel("Aspect Ratio α", color="white", fontsize=12)
ax2.set_ylabel("Trapping Ratio $C_R / C_D$", color="white", fontsize=12)
ax2.set_title("Trapping Ratio vs Aspect Ratio", color="white", fontsize=13, fontweight="bold")
ax2.legend(fontsize=10, facecolor="#1e293b", edgecolor="#334155", labelcolor="white")
ax2.set_facecolor("#0f172a")
ax2.tick_params(colors="gray")
ax2.grid(True, alpha=0.2, color="#334155")
for spine in ax2.spines.values():
    spine.set_color("#334155")

# ─── Plot 3: Gaussian plume with images ───
ax3 = axes[2]
H = 20.0  # building height
z_s = 0.5  # source at 0.5 m (tailpipe)
sigma_z_vals = np.array([2, 5, 10, 15])
z = np.linspace(0, H, 200)

for sig_z in sigma_z_vals:
    C_z = np.zeros_like(z)
    for n in range(-5, 6):
        C_z += np.exp(-(z - z_s - 2*n*H)**2 / (2*sig_z**2))
        C_z += np.exp(-(z + z_s - 2*n*H)**2 / (2*sig_z**2))
    C_z /= C_z.max()
    ax3.plot(C_z, z, linewidth=2, label=f"$\\sigma_z$ = {sig_z} m")

ax3.set_xlabel("Normalised Concentration", color="white", fontsize=12)
ax3.set_ylabel("Height z (m)", color="white", fontsize=12)
ax3.set_title("Vertical Profile\n(Method of Images)", color="white", fontsize=13, fontweight="bold")
ax3.legend(fontsize=9, facecolor="#1e293b", edgecolor="#334155", labelcolor="white")
ax3.set_facecolor("#0f172a")
ax3.tick_params(colors="gray")
ax3.grid(True, alpha=0.2, color="#334155")
for spine in ax3.spines.values():
    spine.set_color("#334155")

fig.patch.set_facecolor("#0f172a")
plt.tight_layout()
plt.savefig("output.png", dpi=120, bbox_inches="tight", facecolor="#0f172a")
plt.close()

# ─── Print summary ───
print("=== OSPM Concentration Analysis ===")
print(f"Canyon width W = {W} m, Background C_B = {C_B} μg/m³")
print(f"Line source q_L = {q_L} g/m/s, Rooftop wind u_H = {u_H} m/s\n")
for a in [0.3, 0.5, 1.0, 1.5, 2.0, 3.0]:
    idx = np.argmin(np.abs(alpha_arr - a))
    print(f"α = {a:.1f}: C_D = {C_direct[idx]:.1f}, C_R = {C_recirc[idx]:.1f}, "
          f"C_total = {C_total[idx]:.1f} μg/m³, trap ratio = {trap_ratio[idx]:.2f}")
    if C_total[idx] > 40:
        print(f"  ⚠ EXCEEDS WHO annual NO₂ limit (40 μg/m³)")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

16.2.6 Key Takeaways

The advection-diffusion equation $\partial C / \partial t + \mathbf{u} \cdot \nabla C = \nabla \cdot (K\nabla C) + S$ governs pollutant transport in street canyons.
The eddy diffusivity $K(z)$ has a parabolic profile vanishing at ground and roof, from mixing length theory.
The method of images enforces reflecting boundary conditions at the canyon walls.
OSPM decomposes concentration as $C = C_D + C_R + C_B$; the trapping ratio $C_R/C_D \sim \alpha$ grows linearly with canyon depth.
Deep canyons with $\alpha > 1$ routinely exceed WHO limits for NO₂ even under moderate traffic loads.

← Canyon Flow Regimes Turbulence & Health Metrics →