Impact of Uncertainties in the Atmospheric Boundary Layer Height on the Numerical Simulation of Chemical Species

Abstract

Uncertainties in the atmospheric boundary layer height (ABL-H) are associated either with the uncertainties in meteorological fields or the definition of the ABL-H. When the ABL-H is involved in ABL parameterization schemes these uncertainties can have impacts on numerical simulation of chemical species. The impacts are examined numerically by employing a 1-D diffusion model with the K-profile scheme and 3-D air quality forecast model with the turbulent kinetic energy (TKE)-based scheme. The dependency on the ABL-H in the two schemes is very different. Sensitivity tests with the 1-D model show that the over-estimated/under-estimated ABL-H leads to the decrease/increase of concentration of tracers under both stable and unstable conditions due to the increase/decrease of the volume of tracers within the ABL. Under unstable conditions, the over-estimated/under-estimated ABL-H also enhances/weakens the vertical diffusivity, and leads to the decrease/increase of the concentration of tracers for the negative/positive vertical gradient of tracer. These impacts result in a little change with height within the ABL. The impacts of the ABL-H through the counter-gradient term and entrainment flux are much smaller than the impact through changing the volume of the ABL and vertical diffusivity in the ABL. Sensitivity tests with the 3-D numerical model with the TKE scheme show that over urban areas the over-estimated/under-estimated ABL-H leads to the increase/decrease of TKE in the whole ABL but does not always lead to stronger/weaker diffusivity and lower/higher concentrations of chemical species.

1 Introduction

The atmospheric boundary layer height (ABL-H) is a parameter describing the vertical extent of turbulent mixing. It is a diagnostic variable in numerical models whose value is determined mostly by the potential temperature profile or Richardson number. The evolution of inversion layers affects the ABL-H, which is what determines the dilution of chemical species within the ABL. Therefore, the ABL-H has a great impact not only on the transport of heat and momentum, but also on the transport and dispersal of chemical species within the ABL.

Culf (1992) investigated the influence of the ABL-H on the averaged CO$_{2}$ concentration over the ABL over the tropical forest. By analysing the results of a simple 1-D bulk model and observations, they found that while daytime CO$_{2}$ mean concentration is controlled almost entirely by mixed layer growth, with the surface flux playing a very small role, nocturnal CO$_{2}$ mean concentrations are much more sensitive to differences in both boundary layer height and surface flux of CO$_{2}$. However, the work of de Arellano et al. (2004) shows that the mixed layer growth is more important in affecting averaged CO$_{2}$ concentration than the surface uptake only in the morning. The surface uptake becomes dominant at other times.

The relationship between the ABL-H and the concentration of air pollutants has been widely investigated. Su et al. (2018) and Su et al. (2020) examined the role of the ABL-H in pollution events. Their results show that there is a negative correlation between the ABL-H and near surface air pollutants. The shallow ABL can lead to severe air pollution episodes even when the anthropogenic emissions were significantly reduced during the lockdown due to the Covid-19 pandemic in 2020. Long term observations suggest the negative correlation between the ABL-H and diffusion capacity of air pollutants (Xiang et al. 2019; Lee et al. 2019). By examining five year observation data sets, Yuval et al. (2020) found that the correlation is different in different types of the ABL: convective boundary layer and ground-based residual layer often lead to negative correlation between the ABL-H and concentrations of NO$_{x}$ and PM$_{2.5}$ and stable boundary layer leads to positive correlation. Schafer et al. (2006) and Wagner and Schafer (2017) examined the correlation over urban areas based on multi-year observations and found that the correlation varies with geographic locations (e.g., urban vs. rural) as well as the type of chemical species.

To understand the relationship, the 1-D bulk model has been employed by Raupach (1991), Denmead et al. (1996), Culf et al. (1997), Dang et al. (2011) and Pino et al. (2012) to examine the impacts of the ABL-H, surface flux, and entrainment on the mean CO$_{2}$ concentration (defined as the vertically integrated mixing ratio from ground to the top of ABL divided by the ABL-H). The results of Pino et al. show analytically the dependency of the mean concentration on the ABL-H under the fluxes from upper and lower boundaries. If the flux from the upper boundary is very small, the mean CO$_{2}$ concentration is inversely proportional to the ABL-H. These results, however, don’t show the impact on concentration at different heights within the ABL. Reen et al. (2014) use a dispersion model to investigate the impact of the ABL-H on the concentrations of air pollutants.

The evolution of the ABL-H is driven by many factors including surface heat flux, weather patterns, advection of meteorological fields, surface roughness and land coverage etc. (Tennekes 1973; Culf 1992). It has long been recognized that there are big uncertainties in the ABL-H which are associated with not only the uncertainties in meteorological variables involved in the definition of the ABL-H, such as temperatures, winds and humidity but also different definitions of ABL-H. The papers of Seidel et al. (2010) and Seidel et al. (2012) show that the uncertainties in meteorological fields can produce more than 50% uncertainty in shallow boundary layers and 20% in deeper boundary layers. In addition, ABL-H with different definitions can have very different values even under the same meteorological conditions (Ren et al. 2020).

Given the large uncertainties in the ABL-H and the relationship between the ABL-H and chemical species, it is important to investigate the possible impacts of the uncertainties on the numerical prediction of concentrations of chemical species if the ABL-H is involved in ABL parameterization schemes. Understanding the impacts can help understand the different concentrations of chemical species simulated with different ABL-H definition or with the same ABL-H definition but different ABL parameterization scheme. More importantly, it can help quantify the uncertainties in air quality prediction associated with the ABL-H uncertainties, and find ways to improve air quality forecast.

Dependency on the ABL-H is different in different ABL parameterization schemes. Among the schemes the K-profile and turbulent kinetic energy (TKE)-based parameterization schemes are two important ones with very different dependency on the ABL-H. In the K-profile scheme the ABL-H involves in both the convective velocity and diffusion coefficient under unstable conditions. In the TKE-based scheme the ABL-H involves in the convective velocity which appears in the lower boundary condition. Therefore the impacts of the ABL-H would be different in the different parameterization scheme. In this work, the two schemes are employed to examine the impacts of the ABL-H uncertainties on numerical simulations of concentrations of chemical species through changing the diffusion coefficient. Unlike previous studies with the 1-D bulk model which consider only the impact on averaged concentration, our work will examine the impact on concentration at different heights.

This paper is organized as follows. A 1-D model with the K-profile parameterization scheme is employed in Sect. 2 to examine the impacts of uncertainties in the ABL-H on vertical diffusivity, tracer’s concentration, counter-gradient term and entrainment flux. In Sect. 3, the role of the ABL-H in the TKE-based parameterization scheme is discussed, and the sensitivities of concentration of chemical species to the uncertainties in the ABL-H are investigated with Global Environmental Multiscale Model-Modelling Air Quality and Chemistry (GEM-MACH). Summary and discussion are given in the last section.

2 Impacts in the K-Profile Parameterization Scheme

The K-profile scheme is one of important ABL parameterization schemes which has been widely employed in numerical models to describe the turbulence contributions to the motions within the ABL. In the scheme, the ABL-H is involved in both convective velocity scale and diffusion coefficient computations under unstable conditions. In this section, a 1-D model with the K-profile scheme is employed to examine the impacts of the uncertainties in the ABL-H on vertical diffusivity, tracer’s concentration and counter-gradient term by comparing the benchmark results with scenario results with the over- and under-estimated ABL-Hs.

2.1 Impact on the Diffusion Coefficient

For unstable conditions (mostly during the daytime), the diffusion coefficient in the K-profile closure scheme proposed by Holtslag and Boville (1993) is computed from the following equation:

$$\begin{aligned} K=kw_{t}z\left( 1-{z\over h}\right) ^2, \end{aligned}$$

(2.1)

where $k=0.4$ is the Karman constant, h is the ABL-H. For stable ABL, $w_{t}=u_{*}/\phi _{T}$, where $\phi _{T}$ is the dimensionless vertical temperature gradient, $u_{*}$ is the friction velocity. For unstable ABL, $w_{t}=w_{m}/Pr$ where Pr is the turbulent Prandtl number,

$$\begin{aligned} w_{m}=(u_{*}^3+c_{1}w_{*}^3)^{1/3}, \quad w_{*}^3={gh\overline{\theta 'w'}|_{s}\over \theta _{vs}}\equiv hB_{0}, \end{aligned}$$

(2.2)

$c_{1}$ is a constant, g is the gravity, $\overline{\theta 'w'}|_{s}$ and $\theta _{vs}$ are heat flux and potential temperature at surface, respectively. Because $w_{t}$ depends on the ABL-H the dependence of K on h is highly nonlinear, and the increase of the ABL-H leads to the increase of K at a given height. In the scheme the ABL-H is involved explicitly in computing the diffusion coefficient under unstable conditions, and diffusivity vanishes at the top of the entrainment zone regardless of the definition of the ABL-H. These two important features don’t appear in the TKE-based scheme.

At midday, $h\sim 1000$ m, $w_{*}\gg u_{*}$. Therefore, at the lower ABL ($h\gg z$), the corresponding change of K to the change of h is:

$$\begin{aligned} {\delta K\over K} \sim {\delta h\over 3h}. \end{aligned}$$

(2.3)

For a 30% change of $\delta h/h$, $\delta K/K$ is about 10%.

In the sensitivity tests, the evolution of the ABL-H is prescribed by:

$$\begin{aligned} h(t)=\left\{ \begin{array}{ll} 100hf+(1000hf-100)\cos \left( {t-12\over 12}\pi \right) &{} 6hr< t < 18hr \\ 100hf &{} 6hr\ge t\ge 0, \quad t\ge 18hr, \end{array} \right. \end{aligned}$$

(2.4)

where hf is a factor for sensitivity test, $hf=1$ for the benchmark case, 1000hf m is the maximum ABL-H at $t=12$. The ABL is assumed to be unstable between 6:00 AM and 6:00 PM. The convective vertical velocity $w_{*}$ is computed using $w_{*}^3=0.015h\cos ((t-12)\pi /12)$, and $w_{m}=(0.3^3+w_{*}^3)^{1/3}$.

The diurnal evolutions of the benchmark diffusion coefficients ($hf=1$) at 30 and 200 m are shown in Fig. 1a, b, respectively. Between 6:00AM to 6:00PM the diffusion coefficient is computed using Eq. (2.1) and is 2 m$^2$ s$^{-1}$ at other times. To examine the impact of uncertainties in the ABL-H on diffusion coefficient, $hf=1.1$, 1.2, 1.3, 1.4 are used to represent the over-estimated ABL-H. The differences between the benchmark diffusion coefficient and the scenario values are shown in Fig. 1c, d, respectively. It can be seen that the over-estimated ABL-H leads to the increase of diffusion coefficient, and the differences are roughly close to the difference estimated by Eq. (2.2). At $z=30$ m, the difference at midday is about 3.25 m$^2$ s$^{-1}$ for $hf=1.4$ in Fig. 1c and is about 3.33 m$^2$ s$^{-1}$ according to Eq. (2.3).

The vertical profiles of the benchmark diffusion coefficient and the departures of the diffusion coefficient with $hf=0.8$ and $hf=1.4$ from the benchmark value at hour 12 and hour 14 are shown in Fig. 2. In the lower ABL, the over-estimated/under-estimated ABL-H leads to the increase/decrease of diffusion coefficient, and the magnitudes of the differences increase monotonically as height increases in the lower ABL. The differences can be comparable with the benchmark diffusion coefficient around 500 m at hour 12 and at hour 14 when $hf=1.4$. They can be 30% of benchmark diffusion coefficient when $hf=0.8$.

Fig. 1

Diurnal variation of diffusion coefficient K with $hf=1$ at 30 m (a) and at 200 m (b). The differences between K with $hf=1.1$, 1.2, 1.3, 1.4 and the benchmark K at 30 and 200 m are shown in (c) and (d), respectively

Fig. 2

Vertical profiles of benchmark diffusion coefficient (K), difference between diffusion coefficient with $\textrm{hf}=1.4$ (DK1.4) and K (dash-dotted line), difference between diffusion coefficient (DK0.8) (thick dashed lines) with $\textrm{hf}=0.8$ at hour 12 (a) and hour 14 (b)

2.2 Impacts on the Tracer’s Concentration

Because the ABL-H involves in the calculation of the vertical diffusion coefficient in Eq. (2.1), the change of the diffusivity associated with the change of the ABL-H shown above should lead to the change of concentration of tracers. To illustrate the impacts on tracer’s concentration ($\psi$), the following 1-D diffusion model (without counter-gradient) term:

$$\begin{aligned} {\partial \psi \over \partial t}={\partial \over \partial z}\left( K{\partial \psi \over \partial z}\right) , \end{aligned}$$

(2.5)

is used. The work of Ren (2019) shows that this simple model can describe the important features of CO$_{2}$ evolution even with a simple time-dependent diffusion coefficient.

Two kinds of surface emission of CO$_{2}$ shown in Fig. 3 are used in the lower boundary condition for Eq. (2.5). The first kind represents the biogenic CO$_{2}$ flux similar to the flux shown in Yi et al. (2000) with upward flux in the early morning and night time due to respiration, and downward flux during the daytime due to photosynthesis of CO$_{2}$ (Fig. 3a). The second kind represents the flux from industrial activities (Fig. 3b) which is always positive. Since we are only interested in the impacts in the lower part of ABL, flux from top of the ABL is not included for simplicity.

Fig. 3

Diurnal variation of the biogenic surface emission (a) and industrial surface emission (b)

The initial condition of $\psi$ (ppmv) below the top of the ABL used in the computations is:

$$\begin{aligned} \psi (z, t=0)=420-0.03\alpha z, \end{aligned}$$

(2.6)

where $\alpha$ is a (positive) constant measuring the vertical gradient of initial condition. 420 ppmv is used as climatological value of CO$_{2}$ in GEM-MACH.

The corresponding change of tracer’s concentration to the change of the diffusion coefficient can be described by the following equation (Ren et al. 2020):

$$\begin{aligned} \delta \psi =-\int _{0}^{t}\int _{0}^{h}{\partial \psi _{b}\over \partial z'}{\partial G\over \partial z'}\Delta Kdz'd\tau , \end{aligned}$$

(2.7)

where $\psi _{b}$ is the benchmark mixing ratio with $hf=1$, $\Delta K$ is the change of diffusion coefficient, G is the Green’s function of Eq. (7) in Ren and Stroud (2020). Equation (2.7) is derived based on a 1-D model but can be generalized in the 3-D case. The work of Ren and Stroud (2020) shows that in the 1-D model with the height-independent diffusion coefficients the gradient of G is negative near the ground. Thus Eq. (2.7) suggests that in the 1-D model the increase/decrease of K would lead to decrease/increase of $\psi$ for negative $\partial \psi _{b}/\partial z$ and increase/decrease of $\psi$ for positive $\partial \psi _{b}/\partial z$. This is important for understanding the following numerical results based on the 1-D diffusion model.

In order to examine the impact of the ABL-H on concentration, $\psi _{b}$ with $hf=1$ is used as the benchmark results. Because the gradient factor ($\alpha$) can modulate the impact, the benchmark mixing ratios for different $\alpha$ are computed with the two kinds of surface emissions. Figures 4 and 5 show the diurnal variation of the departure of the mixing ratio from 420 ppmv at 30 m for the biogenic and industrial surface fluxes, respectively. Due to the diurnal variation of vertical diffusivity, the diurnal variations in both figures show the strong rectifier effect with a big jump of mixing ratio around 6:00 AM and then a quick drop afterwards. The drop lasts longer in Fig. 4 than in Fig. 5 due to the negative surface flux during the daytime. The two figures also show that compared to the mixing ratio with zero gradient, large gradient factor leads to a sharp reduction of mixing ratio after 8:00 AM.

The uncertainties in the ABL-H cause the change not only in the vertical diffusivity but also in the volume occupied by the tracer. Both factors contribute to the sensitivity of mixing ratio shown in the following numerical results. To examine the sensitivity of $\psi$ to the ABL-H, $hf=0.6$, 0.8, 1.3, 1.4 are used to represent four scenarios with under- and over-estimated ABL-H.

Fig. 4

Diurnal variation of the benchmark mixing ratio at 30 m produced with the biogenic surface flux for a gradient factor $\alpha =0$ (a), $\alpha =0.4$ (b), $\alpha =0.6$ (c) and $\alpha =1$ (d). The ABL-H is set to $h=1$ km

Fig. 5

Same as Fig. 4 but with the industrial surface flux

2.2.1 Impact with the Biogenic Surface Flux

Figure 6 shows the diurnal variations of the mixing ratio differences between the four scenario results and the benchmark results at 30 m for different gradient factor with the biogenic surface flux. It can be seen from the figure that in the early morning and nighttime the over-estimated/under-estimated ABL-H leads to decrease/increase of mixing ratio. Because diffusivity does not change with the ABL-H during this period, the change of mixing ratio is entirely due to the change of the volume of tracer in the ABL and is not sensitive to the gradient factor. Comparison between Figs. 4 and 6 suggests that for $hf=0.6$ and $hf=1.4$, the increase and decrease of mixing ratio can be 60 and 20% of the variations of benchmark values, respectively.

During daytime, the impact of the ABL-H on mixing ratio is sensitive to the gradient factor. For $\alpha =0$, over-estimated/under-estimated ABL-H leads to the increase/decrease of mixing ratio. For large $\alpha$, however, the sensitivity changes: the under-estimated/over-estimated ABL-H leads to the reduction/enhancement of mixing ratio after 6:00AM. The comparison between Figs. 4 and 6 shows that the reduction of variations can be 20% of the mixing ratio variations in Fig. 3 for $\alpha =1$ and $hf=1.4$.

Although the vertical gradient of mixing ratio is negative (or zero for $\alpha =0$) in the initial condition, it becomes positive after 6:00 AM due to negative surface flux. Thus the over-estimated/under-estimated ABL-H leads to the increase/decrease of mixing ratio by enhancing/reducing the vertical diffusivity (Eq. (2.7)), and at the same time leads to the decrease/increase of mixing ratio through diluting/concentrating (with the same surface flux). The results suggest that the former effect is dominant for small $\alpha$ (corresponding to large positive vertical gradient of tracer) and the latter effect is dominant for large $\alpha$ as shown in Fig. 6.

The vertical profiles of the impacts of under-estimated ($hf=0.6$) and over-estimated ($hf=1.4$) ABL-H on concentration at hour 12 are shown in Fig. 7. For the over-estimated case, increase of the gradient factor leads to the decrease of concentration and the impacts have very small variation in height. For the under-estimated case, increase of the gradient factor leads to the increase of concentration, and the magnitude of the impacts decays with height except for large gradient factors.

Fig. 6

Differences between mixing ratios with $hf=0.6$, 0.8, 1.3, 1.4 for a gradient factor $\alpha =0$ (a), $\alpha =0.2$ (b), $\alpha =0.6$ (c), and $\alpha =1$ (d) and benchmark mixing ratio at 30 m. The biogenic surface flux is used in the calculations

Fig. 7

Vertical profiles of the differences between the scenario results and benchmark results at hour 12 for $hf=1.4$ (solid lines) and for $hf=0.6$ (dash lines) with the biogenic surface emission

2.2.2 Impact with the Industrial Surface Flux

Because the industrial surface flux is always positive, the vertical gradient of mixing ratio is always negative. Figure 8 shows that the over-estimated/under-estimated ABL-H leads to the decrease/increase of mixing ratio at 30 m all day regardless of the value of $\alpha$. The over-estimated/ under estimated ABL-H leads to the decrease/increase of mixing ratio by enhancing/weakening the vertical mixing and the dilution/concentration effect. The reductions/enhancements can be compariable to the benchmark mixing ratio minus 429 ppmv shown in Fig. 5. Figure 9 shows that the over-estimated ABL-H ($hf=1.4$) leads to the decrease of concentration in the whole ABL, and the magnitude of the decrease does not change much with height. However the under-estimated ABL-H ($hf=0.6$) leads to the increase of concentration and the magnitude of the increase decreases as height increases except for large gradient factor.

The above numerical results show that uncertainties in the ABL-H affect the tracer’s concentration by changing the volume occupied by the tracer in the ABL, and by changing the vertical diffusivity under the unstable condition. This impact depends on both surface flux and the vertical gradient of tracers and has small vertical variation.

Fig. 8

Same as Fig. 6 but with the industrial surface flux

Fig. 9

Same as Fig. 7 but with the industrial surface flux

2.3 Impact on the Counter-Gradient Term

In ABL parameterization schemes, the turbulent flux of turbulent concentration $\psi '$ is parameterized as:

$$\begin{aligned} \overline{w'\psi '}=-K\left( {\partial \psi \over \partial z}-\gamma _{\psi }\right) . \end{aligned}$$

(2.8)

The counter-gradient term $\gamma _{\psi }$ is used to describe the nonlocal transport associated with large-eddy motions in a convective boundary layer which is defined as (e.g., Holtslag and Moeng (1991)):

$$\begin{aligned} \gamma _{\psi }=a{w_{*}\overline{\psi 'w'}|_{s}\over hw_{m}^2}, \end{aligned}$$

(2.9)

where a is a constant between 7 and 10.

If $B_{0}\sim 1.4\times 10^{-2}$ m$^2$ s$^{-3}$ (corresponding to surface heat flux ($\rho C_{p}\overline{\theta 'w'}|_{s}$) of 500 W m$^{-2}$ according to Eq. (2.2)) at midday, $w_{*}\sim 14^{1/3}\sim 2.4$ m s$^{-1}$ for the ABL-H of 1000 m. Note that $w_{*}$ is set to zero for negative (downward) surface heat flux.

Since $w_{*}h=B_{0}^{1/3}h^{4/3}$, $\gamma _{\psi }$ is proportional to $h^{-4/3}$ in very unstable case ($w_{*}\gg u_{*}$). Thus, the corresponding changes of $w^{*}$ and $\gamma _{\psi }$ to the change of h are:

$$\begin{aligned} {\delta w^{*}\over w^{*}}={1\over 3}{\delta h\over h}, \quad {\delta \gamma _{\psi }\over \gamma _{\psi }}=-{4\over 3}{\delta h\over h}. \end{aligned}$$

(2.10)

Numerical calculations based on Eq. (2.9) show that the magnitude of $\gamma _{\psi }$ is about 0.001 to 0.003. Sensitivity testing results show that The uncertainties in $\gamma _{\psi }$ induced by the uncertainties in the ABL-H have a very small impact ($\sim$ 1 ppmv) on the concentration with both biogenic and industrial surface fluxes (not shown).

2.4 Impacts on Concentration Through Entrainment Flux

Hong et al. (2006) proposed a scheme to describe explicitly the impact of entrainment processes on the tracer’s concentration within the ABL. In the scheme the impact is described by entrainment flux in the following equation:

$$\begin{aligned} {\partial \psi \over \partial t}={\partial \over \partial z}\left[ K\left( {\partial \psi \over \partial z}-\gamma _{\psi }\right) -\overline{(w'\psi ')_{h}}\left( {z\over h}\right) ^3\right] . \end{aligned}$$

(2.11)

Here $\overline{(w'\psi ')_{h}}$ is the flux at the inversion layer. It can be also expressed as $\overline{(w'\theta ')_{h}}\Delta \psi _{e}/\Delta \theta _{e}$, where $\overline{(w'\theta ')_{h}}$ is about $-$ 0.15 times the surface flux of buoyancy, $\Delta \psi _{e}$ and $\Delta \theta _{e}$ are the jump of concentration and potential temperature within the entrainment zone, respectively.

Because the ABL-H presents explicitly in Eq. (2.11), It is used to examine the impact of ABL-H uncertainty on the CO$_{2}$ concentration through entrainment flux. We first examine the impact on the concentration by comparing the benchmark results based on Eq. (2.11) (with zero $\gamma _{\psi }$) (with $\Delta \theta _{h}=3$ K, $\Delta \psi =5$ ppmv) and those based on Eq. (2.5). Figure 10a, b show that the differences of the two are not sensitive to the surface emissions, and the difference at the three heights are very similar after 9:00 AM. The difference reaches 1.7 ppmv at 600 m before 4:00 PM, 1.8 ppmv at 300 m at 5:00 PM and 1.8 ppmv at 6:00 PM. Because the concentration above the ABL-H does not evolve with time, these differences don’t change when the ABL-H is below 300 and 600 m. After 6:00 PM the entrainment flux is zero and the difference at 30 m does not change. The impact of the ABL-H uncertainties on the CO$_{2}$ concentration through the entrainment flux can be examined by comparing scenario concentrations and the benchmark concentration. Figure 10c, d show the their differences at 30 m with biogenic surface emission and industrial surface emission, respectively. Comparisons between these differences and the corresponding differences without the entrainment flux effects shown in Figs. 6d and 8d suggest that the impact associated with the entrainment flux uncertainty on concentrations is less than 1.8 ppmv much smaller than the impacts associated with other factors discussed in Sect. 2.2.

Fig. 10

Differences between the benchmark concentrations with and without the entrainment flux effect at 30, 300 and 600 m for biogenic surface emission (a) and industrial surface emission (b). Differences between scenario concentrations and benchmark concentration at 30 m with $\alpha =1$ for biogenic surface emission (c) and industrial surface emission (d)

3 Impacts of the ABL-H in the TKE Parameterization Scheme

3.1 The Role of the ABL-H in the TKE-Based Scheme

The TKE scheme is another ABL parameterization scheme widely used in numerical models (Mailhot and Benoit 1982; Benoit et al. 1989; Belair et al. 1999). The evolution of TKE (E) is described by the following equation (Appendix A):

$$\begin{aligned} {\partial E\over \partial t}=BE^{1/2}-CE^{3/2}+{\partial \over \partial z}\left( K_{M}{\partial E\over \partial z}\right) , \end{aligned}$$

(3.1)

where $K_{M}$ is the diffusion coefficient for momentum, B can be positive or negative depending on the local Richardson number, and C is always positive. The first term on the right represents the combined effect of buoyancy and shear production term. The second term is associated dissipation.

Lower boundary condition for the equation is:

$$\begin{aligned} E(z=0, t)=\left\{ \begin{array}{ll} 3.75u_{*}^2 &{} \quad \textrm{stable} \\ 3.75u_{*}^2+0.2w_{*}^2 &{} \quad \textrm{unstable}. \end{array} \right. \end{aligned}$$

(3.2)

Once TKE is available, the diffusion coefficient for momentum $K_{M}$ is obtained from the following equation:

$$\begin{aligned} K_{M}={c\lambda \sqrt{E}\over \phi _{M}}, \end{aligned}$$

(3.3)

where $c=0.516$, $\lambda =\min [k(z+z_{0}),200]$ is Blackadar’s mixing length, $z_{0}$ is the roughness length, $\phi _{M}$ is the stability function which is a (nonlinear) function of Richardson number.

Unlike the K-profile scheme in which the relationship between K and ABL-H can be described analytically and K is zero at the top of the ABL, the connection between $K_{T}(=K_{M}/Pr)$ and ABL-H in the TKE scheme cannot be explicitly described, and there is no guarantee that $K_{T}$ vanishes at the top of the ABL.

Because the ABL-H is contained in $w_{*}$ it can affect the vertical diffusivity indirectly through the lower boundary condition for the TKE equation. This is very different from the way that the ABL-H affects the vertical diffusivity in the K-profile scheme in which the ABL-H involves in the diffusion coefficient explicitly within the ABL. In addition, because the uncertainties in ABL-H can also affect the stability function in Eq. (3.3), their impact on diffusion coefficient can be very complicated.

Because Eq. (3.1) is a nonlinear equation with time-dependent coefficients, its solution is not available in general. To illustrate the impact of the ABL-H on TKE the solution of the perturbed TKE equation with constant coefficients is derived (Appendix B) from which the change of TKE corresponding to the change of the ABL can be obtained as:

$$\begin{aligned} \delta E={c\lambda \over 2}\int _{0}^{t}{\kappa (\tau )\delta h(\tau )K_{M}\over \sqrt{\pi (t-\tau )^3}}\exp \left[ -\left( \Lambda +{K_{M}^2\over 4(t-\tau )}\right) z\right] d\tau , \end{aligned}$$

(3.4)

where $\Lambda$ is defined in (B.7). Equation (3.4) shows that the TKE increases/decreases for the over-estimated/under-estimated ABL-H and the magnitude of the changes decays with height. However, because $\delta K_{M}$ is proportional to $E^{-1/2}\delta E$ and depends on $\phi _{M}$, the over-estimated/under-estimated ABL-H may not lead to the increase/decrease of $K_{M}$, and $\delta K_{M}$ may not decay with height.

3.2 GEM-MACH

The analyses in Sect. 3.1 are based on the constant coefficient assumption. To examine the impact of the ABL-H on chemical species in more realistic meteorological and chemical environment, numerical models with the TKE-based scheme should be employed. GEM-MACH is a multi-scale chemical weather forecast model developed by Environment and Climate Change Canada. Because GEM-MACH is composed of dynamics, physics, and on-line chemistry and transport modules and uses the TKE parameterization scheme, it is a better modeling system than the dispersion model employed by Reen et al. (2014) to examine the impacts of the ABL-H. In this work, a high resolution (2.5-km) GEM-MACH coupled with the town energy balance model as its urbanization model is employed to examine the impacts of the uncertainties in the ABL height on the concentrations of chemical species. This high resolution GEM-MACH has been evaluated extensively for both meteorological and chemical fields (Leroyer et al. 2011, 2018; Ren et al. 2020) over large urban areas. The surface emissions of chemical species are produced by applying the spatial surrogate fields and temporal profiles to transform the emissions inventories into model-ready emission maps. Details of creating emission data sets for GEM-MACH can be found in Stroud et al. (2020).

In Ren et al. (2020) urban effects on meteorological and chemical fields over urban centers of Toronto, Chicago, New York City and Detroit have been examined. Because urban areas are the major sources of air pollutants and have strong convective velocity due to strong surface heat flux the impacts of the uncertainties in the ABL-H over the four urban centers are investigated in this section. The results can help identify the possible source of uncertainties of air quality forecasts, and the possible common feature of the impacts over urban areas.

3.3 Sensitivity Tests

3.3.1 Sensitivity of TKE and Diffusion Coefficient

The impacts of the uncertainties in the ABL-H on chemical species are examined by the sensitivity test. In this part, three scenario simulations with $h=0.8$H, 1.2H and 1.5H representing the under- and over-estimated ABL-H (H) are carried out. Since the ABL-H affects concentration of chemical species through diffusivity which are associated with the convective velocity in the boundary condition for TKE equation, the impacts of the uncertainties in the ABL-H on the diurnal variations of these variables are investigated by comparing the benchmark GEM-MACH simulation (with H) against the three scenario simulations.

The monthly mean diurnal variations of the corresponding TKE differences averaged over the four urban centers at the lowest model level (22 m) are shown in Fig. 11. It can be seen from the figure that the pattern of the diurnal variation over the four urban centers are similar. over-estimated ABL-H leads to the increase of TKE, and the differences reach a maximum at 2:00 PM when the differences of $w^{*}$ between the scenario and benchmark simulations also reach a maximum. The magnitude of the difference is roughly proportional to the magnitude of the ABL-H difference.

Fig. 11

Diurnal variation of the monthly mean TKE differences at the lowest model level between the three scenario simulations and benchmark simulations averaged over Toronto (a), Chicago (b), New York City (c) and Detroit (d)

In GEM-MACH, the averaged diffusion coefficient over the lowest two model levels (22 and 94 m) is used in computing the concentration of chemical species at the lowest model level. The differences of the averaged coefficient of the two levels between the scenario simulations and benchmark simulation are shown in Fig. 11. The comparison between Figs. 11 and 12 shows that the diffusion coefficient and TKE differences have a very strong positive correlation (Fig. 12).

The vertical profiles of the impacts of the over- and under-estimated ABL-H on TKE and diffusion coefficient over Toronto and New York City are shown in Fig. 13. In the figure, The height of the first model level and the 16th model level are 22 and 30,000 m, respectively. Figure 13a, c suggest that the variations of the impacts on TKE with height over the two areas are similar with positive/negative impact from over-estimated/under-estimated ABL-H agreeing with Eq. (3.4). While the magnitudes of the positive impacts decays with height in general, the magnitudes of the negative impacts increases as height increases in the lower part of ABL over Toronto.

The comparison between Fig. 13b, d shows that the vertical variations of the impacts of the ABL on diffusion coefficient are very different over the two areas. Although the over-estimated ABL-H leads to the increase of diffusion coefficient in the whole ABL over New York City, it leads to the decrease of diffusion coefficient in the upper ABL over Toronto. The magnitude of the impact over New York City is much larger than that over Toronto. The comparison also shows that the under-estimated ABL is different in the two areas. Over Toronto the under-estimated ABL leads to the decrease of diffusion coefficient over the whole ABL over, but it leads to the increase of diffusion coefficient in the lower part of ABL over New York City. According to Eq. (3.3) diffusion coefficient depends not only on TKE but also on the stability function. If the stability function is different the impact on diffusion coefficient and TKE may be different. The very different vertical profiles of the impacts on diffusion coefficient suggest the stability function (in Eq. (3.3)) is also affected by the ABL-H. In addition, because the convective velocity is the product of the ABL-H and surface heat flux, different surface heat flux in the four urban areas can also lead to the different sensitive results.

Fig. 12

Diurnal variation of the monthly mean diffusivity differences averaged over the two lowest model levels (99 and 22 m) between the three scenario simulations and benchmark simulations over Toronto (a), Chicago (b), New York City (c) and Detroit (d). The benchmark diffusion coefficient varies from 2.4 to 70.

Fig. 13

Vertical profiles of TKE difference due to the over-estimated ABL-H ($hf=1.4$, solid line) and the under-estimated ABL-H ($hf=0.8$ dash line) over Toronto (a) and New York City (c) at 2:00 PM, and vertical profiles of diffusion coefficient difference due to the over- and under-estimated ABL-H over Toronto (b) and New York City (d). The height of the lowest model level and the 16th model level are 22 and 30,000 m, respectively. The range of the benchmark TKE is from 0 to 2.5 (m/s)$^2$

3.3.2 Sensitivity of CO Mixing Ratio

The change of diffusion coefficient due to uncertainties in the ABL-H can affect all the chemical species in GEM-MACH. In this part we focus on the impact on CO in July, 2015.

The diurnal evolution pattern of CO emission in summer is similar to that of the industrial CO$_{2}$ emission. The CO emission increases rapidly around 7:00 AM from 40 to 134 (g/s), and continues to increase slowly to 171 (g/s) at 3:00 PM before it starts to decrease. The impact of the ABL-H on the mixing ratio of CO can be seen from Fig. 14 which shows the monthly mean diurnal variation of the mixing ratio differences between the scenario simulations and the benchmark simulation at the lowest model level. The impacts over the four areas are quite different. The magnitudes over Toronto and New York City are much larger than over other two areas. Over both Toronto and New York City the over- and under-estimated ABL-Hs lead to the decrease of CO concentration in the afternoon in general. These results are very different from those with the K-profile scheme. The maximum impact is about 15 $\upmu$g/kg which is about 5% of benchmark concentration. The impacts of the over- and under-estimated ABL-Hs have similar evolution pattern over each area, but have different evolution patterns over the two areas. Over Chicago and Detroit the over-estimated ABL-H also leads to the decrease of CO concentration in general.

Although NO$_2$ is involved chemical reactions during the daytime and the change of the ABL-H affects all the chemical species involved in the chemical reactions, the impact of the ABL-H on NO$_2$ still has a similar evolution pattern to that of CO (not shown). The maximum impact can be 0.8 $\upmu$g/kg over Toronto and New York City which is about 10% of benchmark concentration.

Since the impacts on TKE are large around 2:00 PM, the vertical profiles of the impacts of uncertainties at 2:00 PM are used to show the variations of the impacts within the ABL. Figure 15 shows that the variations are different over different areas. While the over-estimated ABL-H leads to smaller concentration over Toronto, Chicago and Detroit in almost the whole ABL, it leads to larger concentrations in the middle of ABL over New York City. The under-estimated ABL-H leads to smaller concentration in almost all the ABL over Toronto and Detroit, but leads to larger concentration over Chicago in almost all the ABL and in the middle ABL over New York City. Among the four areas the impacts of both over- and under-estimated ABL-H over New York City are much larger than over other area with the ABL. According to Eq. (2.7) these differences are due to the combined effects of different sensitivity of diffusion coefficient and the vertical gradient of benchmark mixing ratio.

Fig. 14

Diurnal variation of the monthly mean differences of the mixing ratio of CO at the lowest model level between the three scenario simulations and benchmark simulations over Toronto (a), Chicago (b), New York City City (c) and Detroit (d)

Fig. 15

Vertical profiles of CO mixing ration differences at 2:00 PM due to over-estimated ABL-H ($hf=1.5$ solid line) and under-estimated ABL-H ($hf=0.8$ dash line) over the urban area of Toronto (a), Chicago, (b) New York City (c), and Detroit (d)

4 Summary and Discussion

Due to uncertainties in meteorological fields and different definitions of the ABL-H there exist large uncertainties in the ABL-H. If the ABL-H is involved in the ABL parameterization scheme employed in numerical models, the uncertainties can affect the concentration of chemical species within the ABL by changing the vertical diffusivity and the volume of tracer in the ABL. These impacts on the simulation of concentration of chemical species are examined by using K-profile and TKE-based parameterization schemes.

In the K-profile scheme, the ABL-H is involved directly in the calculation of diffusion coefficient in a highly nonlinear way. The numerical results based on a 1-D diffusion model with the K-profile scheme show the over-estimated/under-estimated ABL leads to larger/smaller diffusion coefficient and the magnitude of the impacts increases as height increases in the lower part of the ABL but decreases in the upper part of ABL. It also leads to decrease/increase of tracer’s concentration by concentrating/diluting the air within the ABL in the early morning and nighttime. In an unstable ABL, the sensitivity of the concentration of tracers to the ABL-H depends not only on the value of the ABL-H but also on the surface flux and vertical gradient of tracers. For negative surface flux, the over-estimated ABL-H leads to the increase of concentration if the (negative) gradient is weak, but leads to the decrease of concentration if the gradient is strong. However, for positive surface flux, the over-estimated ABL-H always leads to the decrease of concentration for negative gradient of tracers. The variation of the magnitudes of the impacts with height is small in the lower part of the ABL. It is found that the uncertainties in counter-gradient term and the entrainment flux term induced by the ABL-H uncertainty have small impacts on concentration.

In the TKE-based scheme, the ABL-H is involved in the lower boundary condition for the TKE equation under the unstable condition. The sensitivity tests show that while the over-estimated/under-estimated ABL-H leads to large/small TKE in the whole ABL it does not lead to large/small diffusion coefficient in the whole ABL over the four urban areas. Over urban area in Toronto, the over-estimated ABL-H can lead the decrease of diffusion coefficient, and the under-estimated ABL-H can lead to increase of diffusion coefficient over urban area of New York City. The sensitivity results also show that the large/small diffusion coefficient does not always lead to small/large concentration of chemical species during the daytime within the ABL.

The comparison between the sensitivity results with the K-profile and TKE-based schemes show that the uncertainties in the ABL-H have drastically different impacts on both diffusion coefficient and concentration of chemical species. While the over-estimated/under-estimated ABL-H leads to the increase/decrease of diffusion coefficient in the whole ABL in the K-profile scheme, it does not always lead to the increase/decrease of diffusion coefficient in the whole ABL. This difference reflects the different role of the ABL-H in the TKE-based scheme.

With the industrial surface flux, the over-estimated/under-estimated ABL-H leads to the decrease/increase of CO$_{2}$ concentration within the whole ABL in the K-profile scheme. In the TKE scheme, however, the over-estimated ABL-H leads to the decrease of CO concentration only in the few lowest model levels and leads to the increase of concentration at some other heights. The under-estimated can also lead to the decrease of CO concentration. Due to the lack of 3D meteorological fields such as wind and temperature, these are not fair comparisons. Because Green’s functions in the 1-D and 3-D models can be different, and the difference may lead to different sensitivity results according to Eq. (2.7). A better way to compare the sensitivity results is to compare the results of GEM-MACH simulations with the K-profile scheme against the TKE-based scheme. Because the K-profile scheme is not implemented in GEM, such comparison is not available. But it will be out future work.

The work of Su et al. (2020) shows that the low ABL-H can lead to enhancement of air pollutants during the pandemic period. Our results with the K-profile also show the strong anti-correlation in the sensitivity results. Although the monthly mean benchmark results of GEM-MACH also show strong anti-correlation between 8:00: AM and 12:00 PM, the sensitivity results don’t show the much changes in the correlation. This is because unlike the K-profile scheme in which diffusivity is zero at the inversion layer, the change of the ABL-H in the tests has a very small impact on diffusivity at the inversion layer.

Our results show different impacts of the ABL-H uncertainty in the K-profile and TKE-based schemes. These results are not enough to be used as the criteria for assessing the performance of each scheme. This is because that the ABL parameterization scheme contains many parameters such as surface heat flux, stability function, mixing length, ABL height etc. The performance of the scheme depends on the accuracy of all these parameters.

The theoretical results of this work provide qualitative descriptions of the relationship between the uncertainty in the ABL-H and its impact on concentration. It is difficult to use them for qualitatively estimating the ABL-H uncertainty based on the impact or the magnitude of the impact based on the ABL-H uncertainty. Equation (2.7) contains the ABL-H, diffusion coefficient K, vertical gradient of benchmark concentration, and Green’s function. All these variables have strong variations in time and height. In addition, Eq. (2.7) contains integration over time and height. Equation (3.4) is derived based on the constant coefficient assumption.

Although the mechanism of the impact of ABL-H uncertainties is the same in both urban and rural areas, our results have some implications in urban chemical environment studies. Due to the change of thermal and radiative properties by urban surfaces including roads, walls and roofs, the surface heat flux is stronger in urban areas than in rural areas. Because the convective velocity is the product of the ABL-H and surface heat flux, large surface heat flux can amplify the impact of ABL-H uncertainty. Thus, the simulation results of concentration are more sensitive to the uncertainty in the ABL-H in urban areas than in rural areas.

Appendices

Appendix A: TKE Equation

The TKE is defined as $(\overline{u^{'2}+v^{'2}+w^{'2})}/2.$ The evolution of TKE is described by the following equation (Mailhot and Benoit 1982):

$$\begin{aligned} {\partial E\over \partial t }=S+F+T+D, \end{aligned}$$

(A.1)

where:

$$\begin{aligned} S= & {} -\overline{u'w'}{\partial \bar{u}\over \partial z}-\overline{v'w'}{\partial \bar{v}\over \partial z}, \end{aligned}$$

(A.2)

$$\begin{aligned} F= & {} \overline{w'b'}, \end{aligned}$$

(A.3)

$$\begin{aligned} T= & {} {\partial \over \partial z}\left( \overline{w'E'}+{1\over \rho _{0}}\overline{w'p'}\right) , \end{aligned}$$

(A.4)

$$\begin{aligned} D= & {} -\mu \overline{|\nabla \times \textbf{u}|^2}, \end{aligned}$$

(A.5)

are shear production term, buoyancy flux term, transport and pressure work (redistribution term) and dissipation, respectively.

Observations show that S and D are dominate terms at nighttime. While S is positive, D is negative. Both decay rapidly with height and become very small at around 0.6H. During the daytime, S and D are also dominate terms below 0.2H. S decays rapidly with height and becomes positive at 0.3H. F is also positive below 0.8H and decays linearly with height. T is negative below 0.4H and positive above. Its value increases with height. Above 0.4H the magnitude of the four terms is within $-$ 0.5 to 0.5 (Lenschow et al. 1988).

In the TKE scheme $F+S$ is parameterized as $BE^{1/2}$, where B can be positive or negative, D is parameterized as $-CE^{3/2}$, where C is positive, and T is parameterized as $\partial (K_{M}\partial E/\partial z)/\partial z$, where $K_{M}$ is positive.

Appendix B: Solution of the Perturbed TKE Equation

To illustrate the impact of the ABL-H on diffusion coefficient we consider the following perturbed TKE equation due to the change of the ABL-H:

$$\begin{aligned} {\partial \delta E\over \partial t}=\alpha (z,t)\delta E+{\partial \over \partial z}\left( \beta (z,t)\delta E+K_{M}{\partial \delta E\over \partial z}\right) , \end{aligned}$$

(B.1)

where:

$$\begin{aligned} \alpha ={B\over 2\sqrt{E}}-{3\over 2}C\sqrt{E}, \quad \beta ={c\lambda \over 2\phi _{m}\sqrt{E}}{\partial E\over \partial z}. \end{aligned}$$

(B.2)

Lower boundary condition for $\delta E$ is:

$$\begin{aligned} \delta E(z=0, t)=\left\{ \begin{array}{ll} 0 &{} \quad \textrm{stable} \\ {0.4\over 3}{w_{*}^2\over h}\delta h\equiv \kappa \delta h &{} \quad \textrm{unstable}. \end{array} \right. \end{aligned}$$

(B.3)

Analytical solution of Eq. (B.1) cannot be obtained in general. Here we consider the special case in which $\alpha , \beta , K_{M}$ are constants. Applying the Laplacian transform to Eq. (B.1), one has:

$$\begin{aligned} K_{M}{d^2\psi \over dz^2}+\beta {\partial \psi \over \partial z}-(p+\alpha )\psi =0. \end{aligned}$$

(B.4)

By letting $\psi =\phi \exp (-\beta z/2K_{M})$ the above equation can be transformed into the following form:

$$\begin{aligned} {d^2\phi \over dz^2}-{1\over K_{M}}\left( p+\alpha +{\beta ^2\over 4K_{M}}\right) \phi =0. \end{aligned}$$

(B.5)

Using the boundary condition and applying inverse Laplacian transform one obtains:

$$\begin{aligned} \delta E=\int _{0}^{t}{\kappa (\tau )\delta h(\tau )K_{M}\over 2\sqrt{\pi (t-\tau )^3}}\exp \left[ -\left( \Lambda +{K_{M}^2\over 4(t-\tau )}\right) z\right] d\tau , \end{aligned}$$

(B.6)

where:

$$\begin{aligned} \Lambda ={\alpha +0.5\beta \over K_{M}}+{\beta ^2\over 4K_{M}^2}. \end{aligned}$$

(B.7)

In this special case the change of $K_{M}$ corresponding change of the ABL-H is:

$$\begin{aligned} \delta K_{M}={c\lambda \over 4\phi _{m}\sqrt{E}}\int _{0}^{t}{\kappa (\tau )\delta h(\tau )K_{M}\over \sqrt{\pi (t-\tau )^3}}\exp \left[ -\left( \Lambda +{K_{M}^2\over 4(t-\tau )}\right) z\right] d\tau . \end{aligned}$$

(B.8)