4.2.1 Potential Temperature

The change of pressure and temperature as a function of altitude is obtained via the hydrostatic and ideal gas laws (see Chapter 3). Furthermore, a relationship between pressure and temperature may be obtained by assuming that air parcels undergo adiabatic changes. Therefore, the law of Laplace applies:

where γ is the adiabatic index of the gas (here, air), also known as the Laplace coefficient (γ = c_p / c_v, where c_p and c_v are the specific heat capacities at constant pressure and volume, respectively). Applying the ideal gas law (see Chapter 3), the following relationship between pressure and temperature is obtained under adiabatic conditions:

For a diatomic gas (which is the case of the atmosphere, since molecular nitrogen and oxygen account for 99 % of the air), γ = 7/5 = 1.4. Thus:

Therefore:

The potential temperature, θ, is defined as the temperature of an air parcel brought back to sea level at a standard pressure (i.e., P₀ = 1 atm) under adiabatic conditions. Therefore:

TP−0.29=θ P0−0.29(4.5)

Thus:

θ=T(P0P)0.29(4.6)

4.2.2 Neutral Atmosphere

A constant vertical profile of the potential temperature is defined as neutral. Let us consider an air parcel with uniform properties (pressure, temperature, and density). Under neutral conditions, this air parcel will remain in equilibrium with its environment as it moves upward or downward, since its temperature will change adiabatically and will, therefore, remain consistent with the temperature of the surrounding air.

4.2.3 Unstable Atmosphere

The atmosphere is unstable if the vertical profile of the potential temperature is negative, i.e., if the actual temperature decreases faster with height than the adiabatic gradient. An air parcel moving upward becomes warmer than that of its surroundings, because its adiabatic change in temperature is less than the decrease of the surrounding atmospheric temperature. Therefore, according to the ideal gas law, the density of this air parcel becomes less than that of the surrounding air. Thus, this air parcel keeps moving upward. Similarly, an air parcel moving downward undergoes an adiabatic change and its temperature change is less than that of the gradient of the surrounding atmosphere. Therefore, the surrounding temperature becomes greater than that of the air parcel. Thus, this air parcel is denser than the surrounding air and it continues to move downward. Vertical motions in an unstable atmosphere are enhanced because a small perturbation in the vertical direction leads to continuous motion in the direction of the initial perturbation (see Figure 4.1). Such unstable air parcels lead to vertical mixing of air pollution, at least up to an altitude where the atmosphere becomes neutral or stable (see Section 4.2.4). An unstable atmospheric layer is called a mixing layer. It occurs generally within the PBL during the day, because solar radiation warms up the Earth’s surface, but does not affect the temperature of the lower atmosphere (since the air absorbs little ultraviolet radiation within the troposphere, see Chapter 5). As a result, the vertical potential temperature gradient becomes negative.

Figure 4.1. Schematic representation of vertical motions of an air parcel for neutral conditions (left figure), unstable conditions (middle figure), and stable conditions (right figure). The solid arrows indicate the initial perturbation of the air parcel, which undergoes an adiabatic change, and the dashed arrows indicate the following motion of the air parcel, which is due to the difference between the temperature of the air parcel and that of the surrounding air (i.e., density difference).

4.2.4 Stable Atmosphere

The atmosphere is stable if the vertical profile of the potential temperature becomes positive. As an air parcel moves upward, its temperature follows an adiabatic change and it becomes lower than that of the surrounding air. Therefore, the air parcel becomes denser than the surrounding air. As a result, it tends to move back downward. Similarly, an air parcel moving downward undergoes an adiabatic change and its temperature becomes greater than that of the surrounding air. Thus, this air parcel, being less dense than the surrounding air, tends to move back upward. Therefore, a stable atmosphere tends to suppress vertical motions of air parcels. Thus, such an atmospheric layer may become stratified, because there are no, or few, transfers between the different layers (strata) of the atmosphere (see Figure 4.1). A stable atmospheric layer may be present aloft, for example above a mixing layer, or at the surface. This latter case occurs at night when the infrared radiation emitted by the Earth’s surface toward the atmosphere leads to a decrease of the surface temperature and, accordingly, of the surface layer temperature. Therefore, a temperature inversion develops: the gradient of the actual temperature becomes positive and the gradient of the potential temperature becomes also positive. (Note that an atmospheric layer may be stable without a temperature inversion, as long as the actual temperature gradient is less than the adiabatic gradient.) Then, the atmospheric layer near the surface is stable at night. This stable layer will gradually become thicker as the surface temperature keeps decreasing. The PBL layer located above this stable layer is called the residual layer. Pollutants located within this residual layer are then isolated from the surface. However, during daytime, as solar radiation warms the surface, thereby leading to the gradual disappearance of the surface-based stable layer, the previously stable and residual layers merge to form a new mixing layer. Pollutants that were formerly present within the residual layer are mixed vertically and may then interact with the surface. The evolution of the atmospheric layer near the Earth’s surface is depicted schematically in Figure 4.2. The stratosphere shows a positive vertical profile of the potential temperature; therefore, it is a stable atmospheric layer with few vertical air motions.

Figure 4.2. Schematic representation of the daily cycle of the planetary boundary layer.

Source: After Stull (1988).

4.2.5 Effect of Moisture

This description of temperature gradients applies to dry atmospheres. If the atmosphere is moist, its vertical temperature gradient under adiabatic conditions will be less than that of a dry atmosphere. The water content, q (i.e., the specific humidity expressed in g of water per g of moist air), corresponds to a relative humidity (%), which is greater when temperature decreases, because the saturation vapor pressure of water decreases with temperature. In other words, for a given quantity (mass) of water per quantity of air (mass or volume), relative humidity increases when temperature decreases. Accordingly, for a given relative humidity, the water content of the atmosphere decreases when temperature decreases.

Humidity

Relative humidity (%) may be converted to absolute humidity using the method proposed by McRae (1980). According to Dalton’s law, the water vapor concentration expressed in parts per million (ppm; 1 ppm = 1 molecule per million molecules of air) is:

[H2O(ppm)]=104RH(%)Ps,H2O(T)P(4.7)

where T is the ambient temperature, P_s,H2O(T) is the saturation vapor pressure of water, and P is the atmospheric pressure. The saturation vapor pressure of water depends on temperature according to the following polynomial formula:

Ps,H2O(T)=P exp(13.3185ωT−1.9760ωT2−0.6445ωT3−0.1299ωT4)(4.8)

where ω_T = 1 – (373.15/T), where T is expressed in K. For example, at P = 1 atm and T = 25 °C, P_s,H2O (298 K) = 0.031 atm; at T = 5 °C, P_s,H2O (278 K) = 0.0086 atm. At a relative humidity of 50 % and a temperature of 5 °C, there are about 4300 ppm of water vapor.

Converting % to (g m⁻³) is performed as follows:

[H2O(g m−3)]=0.18 RH(%)Ps,H2O(T)R T(4.9)

Thus, a relative humidity of 50 % at 1 atm and 5 °C corresponds to 3.4 g of water vapor per m³ of air. The water vapor mass per volume of air is called the volumetric absolute humidity. The water vapor mass concentration at saturation (RH = 100 %) increases with temperature. For example, at 1 atm, there are 9.4, 17, and 30 g of water vapor per m³ at 10, 20, and 30 °C, respectively.

Absolute humidity may also be defined as the mass of water vapor per mass of air. It can be defined per mass of dry air or per mass of moist air. The absolute humidity per mass of dry air is usually called the mixing ratio (by weight) and the absolute humidity per mass of moist air is usually called the specific humidity. The molar mass of dry air is 29 g mole⁻¹; therefore, according to the ideal gas law, the mass in g of one m³ of dry air is: 29 P_{dry air} /(RT). The mixing ratio is thus related to relative humidity as follows:

[H2O(g/g dry air)]=0.18RH(%)Ps,H2O(T)R TR T29 Pdry air=0.62RH(%)100Ps,H2O(T)Pdry air(4.10)

According to Dalton’s law, the total pressure is the sum of the partial pressures of dry air and water vapor:

P=Pdry air+RH(%)100Ps,H2O(T)(4.11)

Therefore:

[H2O (g/g dry air)]=0.62RH(%)100Ps,H2O(T)(P−RH(%)100Ps,H2O(T))(4.12)

The specific humidity, q, is related to relative humidity as follows:

q=1(1.61PPs,H2O(T)100RH(%))−0.61(4.13)

Given that water vapor is a small fraction of the total mass of an air parcel, the mixing ratio and the specific humidity generally have very close values. For example, at T = 25 °C and P = 1 atm, P_s,H2O = 0.031 atm (see previous equation). For a relative humidity of 100 %, the mixing ratio is calculated to be 0.0198 g H₂O/g dry air and the specific humidity, q, is 0.0195 g H₂O/g moist air, i.e., a difference of less than 2 %.

4.3.1 Wind Direction in the PBL

As discussed in Chapter 3, the wind direction in the free troposphere results from the balance between two forces: the pressure gradient force and the Coriolis force (which is due to the Earth’s rotation). In the PBL, another force comes into play: the frictional force, which represents the effect of the friction of the flow on the Earth’s surface. Therefore, the wind direction at the surface results from the balance among these three forces: the pressure gradient force, the Coriolis force, and the frictional force. The frictional force is opposite to the wind direction. Furthermore, the Coriolis force is proportional to the wind speed. Since the wind speed decreases near the surface (tending toward zero at or near the surface, see Section 4.3.2), the Coriolis force also decreases near the surface. Thus, the influence of the pressure gradient force on the wind direction is greater near the surface. The theoretical solution for the northern hemisphere leads to a wind direction near the surface that is 45 ° to the left of the geostrophic wind direction. Thus, the wind direction evolves from that direction near the surface toward the geostrophic wind direction at the height within the PBL where the frictional force becomes negligible. The change in wind direction as a function of height is called the Ekman spiral. The depth of this Ekman layer, where the change in wind direction occurs, depends on the effect of the surface. It is greater for a situation with complex terrain (urban canopy, forest …) than for a smooth surface (sea, snow, ice …). The equations that govern the vertical profile of the wind speed and direction are described next.

4.3.2 Equations of Atmospheric Flow and Turbulence

The Navier-Stokes equations are used to calculate the atmospheric flow. They are expressed in Cartesian coordinates as follows:

∂u∂t+u∂u∂x+v∂u∂y+w∂u∂z−f v=−1ρa∂P∂x+1ρa (∂τxx∂x+∂τyx∂y+∂τzx∂z)∂v∂t+u∂v∂x+v∂v∂y+w∂v∂z+f u=−1ρa∂P∂y+1ρa (∂τxy∂x+∂τyy∂y+∂τzy∂z)∂w∂t+u∂w∂x+v∂w∂y+w∂w∂z      =−1ρa∂P∂z+ 1ρa (∂τxz∂x+∂τyz∂y+∂τzz∂z)−g(4.16)

where the terms τ_ij represent the shear stresses (kg m⁻¹ s⁻²), which are due to the fluid viscosity. For a newtonian fluid, Newton’s law of viscosity applies. It is written for fluid motion in a single direction as follows:

τyx=μv,a∂u∂y(4.17)

where μ_v,a is the dynamic viscosity of the air (kg m⁻¹ s⁻¹). This equation may be seen as the flux in the y direction of the momentum in the x direction. This law, defined for a fluid moving in a given direction (here, x), has been generalized by Navier, Stokes, and Poisson for a fluid moving in three dimensions. Then, the shear stresses for momentum in the x direction are as follows:

τxx=2μv,a ∂u∂x+(μD−23μv,a) (∂u∂x+∂v∂y+∂w∂z)τyx=μv,a (∂u∂y+∂v∂x)τzx=μv,a (∂u∂z+∂w∂x)(4.18)

The other terms for the y and z directions, τ_xy, τ_yy, τ_zy, τ_xz, τ_yz, and τ_zz, are written similarly. Clearly: τ_xy = τ_yx, τ_yz = τ_zy, and τ_xz = τ_zx. In these equations, μ_v,a is the dynamic viscosity of the fluid and μ_D is the dilatational viscosity of the fluid (μ_D = 0 for a monoatomic gas and μ_D ≈ 0 for air).

Boussinesq introduced several assumptions that apply to the PBL:

1. 
   

   – The kinematic viscosity and the molecular thermal conductivity are constant (i.e., their dependence on pressure and temperature can be neglected).

   
   
2. 
   

   – The heat produced by the viscous forces is negligible compared to solar radiation in the heat transfer equation.

   
   
3. 
   

   – The atmosphere is considered to be incompressible, since it is an open system (except at the Earth’s surface).

   
   
4. 
   

   – The fluctuations of the state variables (pressure, temperature, and density) are small compared to their mean values.

   
   
5. 
   

   – The ratio of the pressure fluctuation and the mean pressure is small compared to the ratio of the temperature fluctuation and the mean temperature and compared to the ratio of the density fluctuation and the mean density.

   
   
6. 
   

   – The fluctuations of the density are important only when the density is multiplied by the gravitational constant; indeed, the Earth’s surface introduces a constraint on the system in the vertical direction.

The last two assumptions are commonly called the Boussinesq approximation or hypothesis: they imply that the change in the air density can be neglected except in the term (ρ_a g dz) of the hydrostatic equation (see Equation 3.2).

Using the Boussinesq hypothesis, the continuity equation for momentum is as follows:

∂u∂x+∂v∂y+∂w∂z=0(4.19)

For u, the friction term may be written as follows (assuming that μ_v,a is spatially uniform):

1ρa(∂τxx∂x+∂τyx∂y+∂τzx∂z)=μv,aρa(∂∂x(2∂u∂x−23(∂u∂x+∂v∂y+∂w∂z))+∂∂y (∂u∂y+∂v∂x)+∂∂z (∂u∂z+∂w∂x))(4.20)

Changing the order of the derivations:

1ρa(∂τxx∂x+∂τyx∂y+∂τzx∂z)=μv,aρa((∂2u∂x2+∂2u∂y2+∂2u∂z2)−23∂∂x(∂u∂x+∂v∂y+∂w∂z)+∂∂x(∂u∂x+∂v∂y+∂w∂z))(4.21)

According to the continuity equation, the last two terms are zero. Since the kinematic viscosity, ν_v,a(m² s⁻¹), is equal to the ratio of the dynamic viscosity and the fluid density, μ_v,a/ρ_a:

1ρa(∂τxx∂x+∂τyx∂y+∂τzx∂z)=νv,a (∂2u∂x2+∂2u∂y2+∂2u∂z2)(4.22)

Thus, the Navier-Stokes equations may be written as follows:

∂u∂t+u∂u∂x+v∂u∂y+w∂u∂z−fv=−1ρa∂P∂x+νv,a ∇2u∂v∂t+u∂v∂x+v∂v∂y+w∂v∂z+fu=−1ρa∂P∂y+νv,a ∇2v∂w∂t+u∂w∂x+v∂w∂y+w∂w∂z   =−1ρa∂P∂z+νv,a ∇2w−g(4.23)

where u, v, and w are the wind velocities (m s⁻¹) in the x, y, and z directions, respectively, f is the Coriolis parameter (s⁻¹), ρ_a is the air density (kg m⁻³), P is the atmospheric pressure (Pa), ν_v,a is the kinematic viscosity of the air (m² s⁻¹), and g is the acceleration due to gravity (m s⁻²).

The Reynolds number, Re, may be used to estimate whether a flow is turbulent or laminar. It represents the ratio of inertial forces and viscous forces, and it is defined as follows:

Re=ulcνv,a(4.24)

where u is the velocity of the flow, l_c is a characteristic length of the flow, and ν_v,a is the kinematic viscosity (ν_v,a = 1.55 × 10⁻⁵ m² s⁻¹ at 1 atm and 25 °C). For a wind velocity of 10 m s⁻¹ and a characteristic length corresponding to the PBL height, say 1,000 m:

A flow is considered to be turbulent when Re > 3,000. Therefore, the atmospheric flow within the PBL is without any doubt turbulent.

In a turbulent flow, the velocity varies randomly around a mean value. Therefore, it may be written as the sum of a mean value and a term that represents the fluctuation around that mean value (according to the statistical representation of Reynolds):

u=u¯+u‘(4.25)

where u¯ represents the mean value of the velocity and u’ represents the random fluctuation around that mean value. The Navier-Stokes equations may then be written as a function of these mean and random terms. Applying the continuity equation for momentum, the Navier-Stokes equations may be rewritten as follows:

∂u∂t+∂u2∂x+∂uv∂y+∂uw∂z−fv=−1ρa∂P∂x+νv,a ∇2u∂v∂t+∂vu∂x+∂v2∂y+∂vw∂z+fu=−1ρa∂P∂y+νv,a ∇2v∂w∂t+∂wu∂x+∂wv∂y+∂w2∂z=−1ρa∂P∂z+νv,a ∇2w−g(4.26)

By decomposing the velocities as the sum of a mean value and a fluctuation around that mean value (Reynolds decomposition), we obtain for u (similar equations may be written for v and w):

∂(u¯+u‘)∂t+∂(u¯+u‘)2∂x+∂(u¯+u‘)(v¯+v‘)∂y+∂(u¯+u‘)(w¯+w‘)∂z−f(v¯+v‘)=−1ρa∂P∂x+νv,a ∇2(u¯+u‘)(4.27)

Averaging over time leads to the Navier-Stokes equation averaged according to the Reynolds decomposition, i.e., RANS for “Reynolds-averaged Navier-Stokes.” Since the mean of the fluctuation term around the mean is by definition zero, the following relationships apply:

u¯¯=u¯u¯‘=0(4.28)

The RANS equation for u is then as follows:

∂u¯∂t+∂u¯2∂x+∂u‘2¯∂x+∂u¯v¯∂y+∂u‘v‘¯∂y+∂u¯w¯∂z+∂u‘w‘¯∂z−f v¯=−1ρa∂P∂x+νv,a ∇2u¯(4.29)

The terms that include the random variables u’, v’, and w’ represent turbulence, and their products are called the Reynolds stresses:

Rxx=u‘2¯ ;  Rxy=u‘v‘¯ ;  Rxz=u‘w‘¯(4.30)

Similar terms may be obtained for v and w. Note that these Reynolds stresses, as defined here, have units of m² s⁻² and not the units of a standard stress, which are kg m⁻¹ s⁻²; the unit of a standard stress would be obtained if the terms on the right-hand sides were multiplied by the fluid density.

There is a set of nine Reynolds stresses, which are called the Reynolds tensor. Because of symmetry between R_xy and R_yx, R_xz and R_zx, and R_yz and R_zy, there are six unknowns in the Navier-Stokes equations, which result from the turbulent nature of the atmospheric flow. New equations must, therefore, be introduced to “close” this system. Several approaches have been proposed to close the RANS equations with an explicit representation of the turbulent terms.

4.3.3 Closure of the Equations of Turbulent Flows

The Boussinesq hypothesis allows one to introduce a relationship between the momentum turbulent term and the mean value of the flow velocity. One introduces the turbulent kinematic viscosity, K_M (m² s⁻¹), which is defined as follows (here for the x and z directions) by analogy with the shear stress (τzx=μv,a∂u∂z) :

−u‘w‘¯=KM∂u¯∂z(4.31)

However, note that this analogy is approximate, because the dynamic viscosity, μ_v,a, is a property of the fluid (here, the air), whereas the turbulent kinematic viscosity is a characteristic of the flow. K_M must be defined. Prandtl introduced the concept of a characteristic mixing length, l_m, which characterizes the length over which a turbulent eddy will retain its properties before mixing with the surrounding fluid. A dimensional analysis implies that the turbulent kinematic viscosity may be expressed as the product of a characteristic velocity, u_* (called the friction velocity), and this characteristic mixing length, l_m:

The friction velocity is defined as the square root of the product of the fluctuations of the horizontal and vertical velocities:

−u‘w‘¯=u*2 (4.33)

The friction velocity represents the amplitude of the vertical turbulent flux of the horizontal momentum within the surface layer. Equations 4.31 to 4.33 lead to the following results:

u*2=u*lm∂u¯∂z; i.e.,  u*=lm∂u¯∂z(4.34)

The mixing length is not constant, but depends on the distance from the wall (here, the Earth’s surface). This relationship uses the von Kármán constant, κ, taken here to be equal to 0.4:

Therefore, the turbulent kinematic viscosity is expressed as follows:

The friction velocity, u* , is on the order of 0.1 to 0.5 m s⁻¹ for standard atmospheric conditions.

4.3.4 Kinetic Energy of the Atmospheric Flow

The specific kinetic energy of the flow (i.e., kinetic energy per unit mass; m² s⁻²) is defined as follows:

Ek=12(u2+v2+w2)=12((u¯+u‘)2+(v¯+v‘)2+(w¯+w‘)2)(4.37)

The mean specific kinetic energy, E_kEk , is defined as follows:

Ek=12(u¯2+v¯2+w¯2)+12(u‘2¯+v‘2¯+w‘2¯)kMKE=12(u¯2+v¯2+w¯2)kTKE=12(u‘2¯+v‘2¯+w‘2¯)Ek=kMKE+kTKE(4.38)

where k_MKE is the specific kinetic energy of the mean flow and k_TKE is the mean specific kinetic energy of the turbulent component of the flow (called TKE for ”turbulent kinetic energy”). An equation may be derived for the change in TKE. This equation contains a term, ε, which represents the TKE mean dissipation rate (m² s⁻³). Several computational fluid dynamics (CFD) models use the TKE equation to close the RANS equations. Among the models that are the most widely used, one may mention the k-l_m, k-ε, and k-ω models (ω is the TKE specific mean dissipation rate). For a more in-depth description of turbulent flow modeling, see, for example, Nieuwstadt et al. (2016).

4.3.5 Vertical Profile of Wind Speed in the Surface Layer

Wind speed changes with height. The vertical profile of wind speed follows a logarithmic law under adiabatic conditions, which results from the assumptions introduced with the first-order closure of turbulence and the use of the Prandtl mixing length (according to Equations 4.34 and 4.35):

∂u∂z=u*lm=u*κ z(4.39)

In the surface layer, the wind speed in a given horizontal direction is given by a differential equation (assuming here for simplicity that the wind direction is in the x direction; i.e., v = 0):

du(z)=u*κ dzz(4.40)

For a viscous fluid, the no-slip assumption is generally applied, i.e., the fluid velocity is zero at the surface. However, such a boundary condition does not apply here since the logarithm of zero is not defined and one must introduce a height at which the wind velocity is zero. This height, z₀, is called the roughness length. Thus, the wind velocity follows the following profile:

u(z)=u*κ ln(zz0)  for z≥z0(4.41)

The height at which the wind speed becomes zero can be estimated experimentally by extrapolating the wind speed in the vertical direction. Thus, the roughness length can be estimated. This roughness length represents the effect of terrain and obstacles (buildings, vegetation, etc.) on the wind speed vertical profile. As a first approximation, it can be estimated to be in the range of 1/30 to 1/10 of the height of the terrain elements. For example, over grassland with grass about 3 cm high, the roughness length would be in the range of 1 to 3 mm; in an urban area with buildings about 10 m high, the roughness length would be in the range of 30 cm to 1 m.

This logarithmic representation of the vertical wind speed profile is idealized, because no momentum is taken into account below z₀. Furthermore, in the case where obstacles affect the atmospheric flow significantly (for example, in a forest or in an urban canopy), a displacement height, d_b, is introduced, which corresponds to a vertical translation of the surface toward the top of the canopy. For example, for a built area, one may define the displacement height empirically as a function of the fraction of the area that is built, λ_b, and of the mean building height, h_b (MacDonald et al., 1998):

db=hb[1−1−λb4λb] (4.42)

For an urban area with a built fraction of, say, 50 % (λ_b = 0.5), d_b = 0.75 h_b. In the limit, d_b tends toward 0 as λ_b tends toward 0 and d_b tends toward h_b as λ_b tends toward 1.

The vertical profile of the wind speed is then defined above a height (d_b + z₀), at which the wind speed becomes zero. If one wants to represent the flow below that height, one must introduce another model of the atmospheric flow within the canopy. For example, in an urban canopy, experimental data have led to the assumption that the vertical profile of the wind speed may be represented by an exponential function (MacDonald, 2000):

u(z)=ξb u(hb) exp(βb(zhb−1))(4.43)

where ξ_b and β_b are empirical parameters, which depend on the built area configuration, u(h_b) is the wind speed at the top of the canopy (which can be taken to be the wind speed above the canopy multiplied by the cosine of the angle of the wind direction and the street-canyon axis), and h_b is the average building height. For example, the following parameterizations may be used (Masson, 2000):

βb=hb2Ws

ξb=1 for hbWs≤13ξb=1+3(2π−1) (hbWs−13) for 13≤hbWs≤23 ξb=2πfor 23≤hbWs (4.44)

where W_s is the mean street width. The term h_b / W_s is sometimes called the aspect ratio of the street canyon (Landsberg, 1981). Near the street surface, the no-slip condition cannot be satisfied by an exponential profile and a logarithmic profile is used. The vertical profile of the wind speed in an urban area may, therefore, be seen as an ensemble of three connected profiles, which are, starting at the surface:

1. 
   

   – a logarithmic profile near the surface

   
   
2. 
   

   – an exponential profile within the urban canopy

   
   
3. 
   

   – a logarithmic profile above the urban canopy

These three profiles must connect at heights that are defined by the canopy configuration and atmospheric flow in order to obtain a continuous vertical wind profile (e.g., MacDonald, 2000).

For non-adiabatic conditions, vertical motions are affected by the temperature vertical profile. For neutral adiabatic conditions, one may define the non-dimensional wind shear number, Φ_M, as follows:

dudzκ zu*=ΦM(4.45)

where by definition, Φ_M = 1, according to Equation 4.39. This approach, which involves the use of dimensionless variables (sometimes called “numbers”), has been widely used to define the characteristics of the PBL. Once the standard variables (vertical fluxes, etc.) have been represented using dimensionless functions, the results are applicable regardless of the absolute dimensions. This approach is, therefore, called “similarity theory.” For the PBL, one refers to the similarity theory of Monin-Obukhov (after the names of the two Russian scientists, Andrei Monin and Alexander Obukhov, who proposed it in 1954). This approach has been useful to estimate non-dimensional variables using experimental data for a wide range of atmospheric conditions. Thus, Businger and co-workers (1971) and Dyer (1974) have independently developed corrected vertical profiles for unstable and stable conditions. These profiles use a dimensionless height, ζ. For stable conditions:

For unstable conditions:

These functions were determined experimentally with the assumption that the von Kármán constant was κ = 0.35. For a value κ = 0.4, the following relationships should be used (Jacobson, 2005):

Stable conditions: ΦM=1+6ζNeutral conditions: ΦM=1Unstable conditions: ΦM=(1−19.3ζ)−1/4(4.48)

In these equations, the dimensionless height, ζ, is defined as (z/L_MO), where L_MO is the Monin-Obukhov length (as discussed in Section 4.4.3). During daytime, the absolute value of L_MO corresponds to the height above which turbulence is generated more by buoyancy (i.e., thermal processes) than by wind shear (i.e., mechanical processes). Its formulation is provided in Equation 4.59. L_MO is negative for unstable conditions and is positive for stable conditions. For neutral conditions, L_MO is infinite and, therefore, ζ = 0.

The vertical profiles of the wind speed are then obtained by integration:

∫0udu=∫z0zΦM(z)u*κ dzz(4.49)

Analytical solutions may be obtained. For neutral conditions, the solution is straightforward and was provided in Equation 4.41. For stable conditions:

u(z)=u*κ[ln(zz0)+(ΦM−ΦM0)](4.50)

where Φ_M0 corresponds to the value of Φ_M at z = z₀. For unstable conditions:

u(z)=u*κ[ln(1−ΦM1+ΦM)−ln((1−ΦM0) (1+ΦM0))+2(arctan(ΦM−1)−arctan(ΦM0−1))](4.51)

For neutral conditions, ζ= 0 and Φ_M = 1; therefore, the two profiles represented by Equations 4.50 and 4.51 then become equivalent to the logarithmic profile given by Equation 4.41.

Figure 4.3 depicts vertical profiles of the horizontal wind speed in a surface layer of 100 m depth for stable conditions (Equation 4.50), neutral conditions (Equation 4.41), and unstable conditions (Equation 4.51). The assumptions regarding the roughness length, the friction velocity, and the Monin-Obukhov length are provided in the figure caption. The vertical gradient of the wind speed is much more pronounced when the atmosphere is stratified, and it becomes rapidly negligible above the surface when the conditions are unstable.

Figure 4.3. Vertical profiles of the horizontal wind speed in the surface layer for stable, neutral, and unstable conditions above 1 m. The assumptions are a roughness length (z₀) of 0.1 m, a friction velocity (u_*) of 0.3, 0.4, and 0.5 m s⁻¹, respectively, and a Monin-Obukhov length (L_MO) of 50 m, infinite, and −20 m, respectively.

4.4.1 Heat Transfer Equation

The heat transfer equation may be written as follows, using the virtual potential temperature, θ_v:

∂θv∂t+u∂θv∂x+v∂θv∂y+w∂θv∂z=θvcp, dry TvdQdt(4.52)

where T_v is the virtual temperature (K), Q is the amount of heat exchanged by the system (J kg⁻¹), c_p,dry is the specific heat capacity of dry air (J kg⁻¹ K⁻¹), and the term on the right-hand side represents the sources and sinks of heat, i.e., the rates of heat transfer due to changes in water phases (condensation, evaporation, freezing, and melting), solar radiation (absorption of solar radiation in the ultraviolet, visible, and infrared ranges and infrared radiation emission by the Earth’s surface and by greenhouse gases), and human activities (heating, transportation, industrial activities, etc.). This equation may be rewritten as follows:

∂θv∂t+∂uθv∂x+∂vθv∂y+∂wθv∂z−θv(∂u∂x+∂v∂y+∂w∂z)=θvcp, dry TvdQdt(4.53)

The Boussinesq assumption implies (see Equation 4.19) that the fifth term is zero.

4.4.2 Sensible Heat Flux

Turbulence affects the atmospheric flow, as well as any associated variable. Thus, temperature is a turbulent variable, which can be decomposed following the Reynolds approach using a mean value and a random component around that mean value. The heat transfer equation in its Reynolds average form is written as follows:

θv=θv¯+θv‘∂θv¯∂t+u¯∂θv¯∂x+v¯∂θv¯∂y+w¯∂θv¯∂z+∂u‘θv‘¯∂x+∂v‘θv‘¯∂y+∂w‘θv‘¯∂z=θv¯cp, dry TvdQdt(4.54)

A first-order closure may be performed using the turbulent heat transfer coefficient, K_H (m² s⁻¹):

∂u‘θv‘¯∂x=−∂∂x(KH∂θv¯∂x) ∂v‘θv‘¯∂y=−∂∂y(KH∂θv¯∂y)∂w‘θv‘¯∂z=−∂∂z(KH∂θv¯∂z)(4.55)

The last term represents heat transfer via the vertical turbulent flux, which transfers heat from the Earth’s surface toward the atmosphere (see Chapter 5). It is related to the sensible heat flux, which is defined as follows (J m⁻² s⁻¹):

FQs=ρa cp, dry  w‘θv‘¯(4.56)

Then, the heat transfer equation is written as follows (using for simplicity the notation θ_vθv for θ¯v and u, v, and w for u¯ , v¯ , and w¯ ):

∂θv∂t+u∂θv∂x+v∂θv∂y+w∂θv∂z−∇KH∇θv=θvcp, dry TvdQdt(4.57)

4.4.3 Variables Characterizing the PBL Stability

In the surface layer, a temperature scale, θ_*, may be introduced, which is related to the sensible heat flux at the surface and which characterizes the production of turbulence via heat transfer:

θ*=−(w‘θ‘¯)0u* (4.58)

Using the friction velocity, which characterizes the production of turbulence via mechanical processes (flow affected by terrain and/or obstacles), and the temperature scale, a length scale, L_MO, may be defined, which is called the Monin-Obukhov length:

LMO=u*2 θvκ g θ*(4.59)

Thus, the Monin-Obukhov length (Monin and Obukhov, 1954) represents during daytime the height above which the production of turbulence via buoyancy dominates the production of turbulence via mechanical processes (i.e., wind shear). A negative value of L_MO implies that the atmosphere is unstable because of heat transfer. A positive value of L_MO corresponds to a stable atmosphere. Under neutral conditions, L_MO tends to infinity, since the heat flux tends toward zero and, therefore, θ_* tends toward zero. The following ranges of L_MO are typical of atmospheric conditions: very unstable conditions when L_MO ranges between 0 and –10 m, unstable conditions when L_MO is less than –10 m, neutral conditions when L_MO has very high values (>100 m in absolute value), stable conditions when L_MO is greater than 30 m, and very stable conditions when L_MO ranges between 0 and 30 m. However, these ranges are rough estimates.

The non-dimensional Richardson numbers are also used to characterize the atmospheric stability (they are named after the English mathematician and physicist Lewis Fry Richardson). The flux Richardson number is expressed as follows:

Rif=gTvw‘θ‘¯u‘w‘¯∂u∂z(4.60)

The gradient Richardson number is expressed as follows:

Rig=gTv∂θ∂z(∂u∂z)2(4.61)

where T_v is in K. The bulk Richardson number, Ri_b, is derived from Ri_g by using discretized vertical gradients of temperature and horizontal wind speed, obtained from meteorological measurements at distinct heights:

Rib=gTvΔθΔz(ΔuΔz)2(4.62)

A positive Richardson number corresponds to a stable atmosphere and a negative Richardson number corresponds to an unstable atmosphere. For a neutral atmosphere, Ri = 0. The Richardson numbers depend on the height at which they are estimated, since the vertical gradients of wind speed and temperature are functions of z. At a given height (or within a given atmospheric layer in the case of Ri_g and Ri_b), they represent the ratio of the heat-generated turbulence (represented by the sensible heat flux or the vertical temperature gradient) and mechanically generated turbulence (represented by the friction velocity and/or the vertical gradient of the horizontal wind speed). Therefore, the Richardson number is qualitatively the inverse of the Monin-Obukhov length (in a dimensionless form) and there is a theoretical relationship between Ri_g (or Ri_b) and the dimensionless form of the Monin-Obukhov length:

zLMO=φ Rig(4.63)

where z is the height at which Ri_g (or Ri_b) is estimated and φ is a dimensionless parameter, which depends on the meteorological conditions; φ is a function of (z / L_MO) in the relationship proposed by Monin and Obukhov (1954). For example, according to experiments conducted over the sea, φ may vary between 4 and 21 depending on the temperature gradient and the wind speed (Hsu, 1989). Other relationships between L_MO and Ri_b have been developed based on theoretical considerations and experimental data, in order to estimate L_MO from the same measurements as those used to estimate Ri_b.

The convective velocity, w_*, is sometimes used to quantify the intensity of vertical motions within the PBL during unstable conditions:

w*=(g zi(w‘θv‘¯)0θv)1/3(4.64)

where g is the acceleration of gravity, z_i is the PBL height, θ_v is the virtual potential temperature, and (w‘θv‘¯)0 represents the sensible heat source at the surface. The convective velocity is on the order of 1 to 2 m s⁻¹ for daytime conditions with strong solar radiation.

For calm atmospheric conditions (u ≈ 0), the friction velocity tends toward zero and, therefore, the Monin-Obukhov length also tends toward zero. Therefore, it is not appropriate as a measure of atmospheric stability if heat-generated convection is important. Under such conditions, one may use the convective temperature scale, θ_*C, which is defined with respect to w_* (instead of θ_*, which is defined with respect to u_*, as shown in Equation 4.58), to characterize the importance of the local thermal convection:

θ*C=−(w‘θ‘¯)0w*(4.65)

4.4.4 Vertical Temperature Profiles

The vertical heat flux can be expressed in a dimensionless manner, introducing the variable Φ_H, similarly to the vertical momentum flux:

dθvdzκ zθ*=ΦH(4.66)

where θ_v is the virtual potential temperature and θ_* is the potential temperature scale defined previously. The virtual potential temperature profile is obtained by integrating in the vertical direction:

∫θv0θvdθv=∫z0zΦH(z)θ*κ dzz (4.67)

For a neutral atmosphere (i.e., under adiabatic conditions), Φ_H = Pr_t, where Pr_t is the turbulent Prandtl number. This dimensionless number is defined as the ratio of the turbulent diffusion coefficients for momentum and heat. Businger and coworkers estimated Pr_t = 0.74 assuming that the von Kármán constant was κ = 0.35. Assuming that κ = 0.4 (i.e., the value most commonly used), then using Pr_t = 0.95 leads to the same results as those obtained by Businger et al. Similarly to the approach used to obtain the vertical profiles of the horizontal wind speed, empirical values have been developed for Φ_H to obtain vertical profiles of the virtual potential temperature under stable and unstable conditions. These profiles are expressed as a function of a dimensionless height, ζ = z/L_MO. The functions of Φ_H developed by Businger and coworkers used the assumption that κ = 0.35. The functions presented next correspond to κ = 0.4 (Jacobson, 2005).

Under stable conditions:

Under unstable conditions:

Then, the vertical profiles of the virtual potential temperature are as follows. Under neutral conditions:

θv(z)=θv(z0)+θ*κ 0.95 ln(zz0)(4.70)

Since there is no sensible heat flux: w‘θ‘¯=0 , therefore, θ_* = 0 and θ_v(z) = θ_v(z₀). Thus, the potential temperature is constant with height under neutral conditions.

Under stable conditions:

θv(z)=θv(z0)+θ*κ [0.95 ln(zz0)+7.8LMO(z−z0)](4.71)

Under unstable conditions:

θv(z)=θv(z0)+θ*κ0.95 [ln[0.95−ΦH0.95+ΦH]−ln[0.95−ΦH00.95+ΦH0]](4.72)

where Φ_H0 corresponds to the value of Φ_H at z = z₀. The vertical profiles of the virtual potential temperature given by Equations 4.71 and 4.72 tend toward a constant value when the atmospheric conditions become neutral, since the Monin-Obukhov length then tends toward infinity and θ_* tends toward zero.

4.4.5 Latent Heat Flux

The latent heat flux corresponds to the water vapor flux associated with evapotranspiration of vegetation and evaporation of water from the surface. This vertical transport of water vapor from the Earth’s surface to the atmosphere results from the heat added to the surface by solar radiation. When the water vapor condenses later in the atmosphere (leading to cloud formation), the condensation process leads to heat production. Therefore, the water vapor flux from the surface to the atmosphere corresponds to a hidden heat flux (the Latin “latens” means hidden).