IV. Planetary Boundary Layer

A. Vertical diffusion - Mellor-Yamada closure scheme

The effects of mixing of heat, momentum and moisture by small scale turbulence is represented by vertical diffusion in the COLA GCM. The mixing coefficients are calculated according to the "level 2.0" closure scheme of Mellor and Yamada (1982). This method assumes a local balance between production and dissipation of turbulent kinetic energy. There are no explicit prognostic variables to describe the planetary boundary layer (PBL); instead, the entire atmosphere is represented in discrete layers which may or may not be part of the PBL. The prognostic equations for atmospheric temperatures and moisture are then coupled to the SSIB equations for the ground surface and canopy, and the system of coupled equations is solved simultaneously with vertical diffusion of heat, moisture, and momentum as given by the Mellor and Yamada (1982) scheme.

In the turbulence closure method, each prognostic variable is expressed as a sum of a large scale (resolved) part and a turbulent (sub-grid) scale part. The vertical fluxes that are expressed as quadratic terms in the turbulent quantities are assumed to be represented by vertical diffusion down the gradient of the large scale quantities, e.g.

where u and w are the turbulent zonal and vertical velocity components, U is the large scale zonal velocity components and K_M is a momentum diffusion coefficient. The overbar denotes averaging on the large (grid) scale, and z represents the vertical direction coordinate. The turbulence closure method involves two assumptions. (1) The diffusion coefficients are expressed as:

where K_M and K_H are the diffusion coefficients for momentum and heat, respectively, l is the master turbulence length scale, q² is the turbulent kinetic energy (so q is the magnitude of the turbulent wind velocity), and S_M and S_H are momentum flux and heat flux stability parameters, respectively. The ratio of S_H to S_M is equal to the ratio of the turbulent flux Richardson number to the bulk (large scale) Richardson number. K_M and K_H are constrained, viz.

and (2) The turbulent kinetic energy (TKE) equation is closed by assuming that production and dissipation of TKE is balanced instantaneously at each point. This assumption leads to a system of equations for l, q, S_M, and S_H that may be solved simultaneously to give vertical diffusion coefficients used to couple the discrete layers of the GCM and the surface layers in SSiB.

B. Gravity wave drag

The parameterization of gravity wave drag (GWD) is based on simplified theoretical concepts and observational evidence. The parameterization consists of determines the drag due to gravity waves at the surface and its vertical variation. We have adopted the scheme of Alpert et al. (1988) as implemented by Kirtman et al. (1992). The parameterization of the drag at the surface is based on the formulation described by Pierrehumbert (1987) and Helfand et al. (1987).

 The momentum flux due to the gravity waves averaged over a grid box is written as

where r is the density of air and V’ and w’ are the fluctuations in the horizontal wind vector and vertical velocity due to gravity waves, respectively. The angle brackets indicate an average over the model grid box. The GWD due to this momentum flux is then an additional body force on the atmosphere. That is,

Pierrehumbert (1987) argued that convective instability and wave breaking at the Earth’s surface are so prevalent that nonlinear effects become very common and cannot be neglected. Based on scaling arguments and the results from numerical experiments, Pierrehumbert (1987) obtained a formula for the surface wind stress due to gravity waves. He concluded that

where F_r = Nh/U is the Froude number, N is the Brunt-Vaisala frequency, U is the surface wind speed, h is the amplitude of the orographic perturbation, l* is the wave length of the monochromatic wave in the direction of the surface wind and t _gw is in the direction of the surface wind.

The equation for the surface wind stress is valid for a wide range of values of the Froude number, in particular for F_r < 0.8 for which the linear theory is valid, and for F_r > 0.8 where nonlinear effects become important. On the other hand, the momentum flux due to the gravity waves above the wave breaking level increases linearly with Froude number even for F_r > 0.8. Assuming that the momentum flux due to the gravity waves does not change with height below the breaking level, the difference between the surface drag, which equals t _gw below the breaking level, and the momentum flux due to the gravity waves above the breaking level is proportional to the wave drag in the intervening layer. Consequently, in the nonlinear range, the momentum deposition within the breaking layer is much larger than in the linear range.

To account for the nonlinear amplification of the surface drag, a base layer is defined to compute t _gw at the surface. The base layer is approximately the lowest one third of the model atmosphere. The equation for the surface stress is computed by using mass weighted variables over the depth of the base layer. In the COLA GCM, l* = 125 km and h equals the subgrid scale standard deviation of the orography. The standard deviation is computed from the Navy 10’ by 10’ orography data by taking departures from the mean orography of the model grid. As in other studies, the computed h is constrained to 400 m for numerical stability. This constraint is based on the assumption that the Coriolis acceleration is at least one order of magnitude larger than the acceleration due to gravity wave drag (see Helfand et al., 1987).

As shown by Eliassen and Palm (1961), t _gw is independent of height in the absence of turbulent dissipation or transience. In such a case, the parameterization deposits the momentum in the top layer of the model and hence there is no body force on the atmosphere in the remaining lower layers because ¶ t _gw/¶ s = 0.

On the other hand, the value of t _gw changes in the vertical when the local value of the wave modified Richardson number, Ri_m = Ri (1-F_r) is less than or equal to 0.25. Wave breaking occurs at levels where amplification of the wave causes the local Froude number to exceed a critical Froude number, F_rc = (1 – 1/4Ri). If such a condition holds, the wave saturation hypothesis of Lindzen (1981) is invoked and wave induced turbulence reduces the magnitude of t _gw. Finally, t _gw vanishes when a critical level, convective instability, or shear instability is encountered. Critical levels occur where V · V_surface < 0. Shear instability, where Ri < ¼, creates turbulence and, according to Eliassen and Palm (1961), the momentum flux vanishes at that level. In order to allow the gravity wave effects to propagate above the base layer, it is assumed