Divergence damping

The divergence damping terms are essentially used to control computational modes (for a hydrostatic model they are used to damp gravity waves). When these terms are not used, noise may develop in the divergence fields (which may eventually contaminate the forecast). They take the form of a diffusion term, and are expressed in general form as

$\displaystyle D_{d\,u}^o$ $\displaystyle = \wp_{o} {\frac{\kappa_{do}}{\mu_d}}
{\frac{\partial^{o-1}}{\partial x^{o-1}}} \left(\nabla\cdot{\bf V}
\mu_d\right)$ (82)
$\displaystyle D_{d\,v}^o$ $\displaystyle = \wp_{o} {\frac{\kappa_{do}}{\mu_d}}
{\frac{\partial^{o-1}}{\partial y^{o-1}}} \left(\nabla\cdot{\bf V}
\mu_d\right)$ (83)

where the damping coefficient, $\kappa_{do}$ (m$^2$ s$^{-1}$), is defined as

$\displaystyle \kappa_{do} = {\frac{k_{d}\, \Delta x^{o} }{2^{o} \, \Delta t}}
\hskip1.in
\left(k_{d} \leq 1\right)$ (84)

where $\Delta x$ represents the grid spacing (m), and $\Delta t$ is the time step (s). Note that $k_{d} \leq 1$ to ensure numerical stability using an explicit time stepping procedure (used in the current version of the model). The leading coefficient simply controls the sign and is defined as

$\displaystyle \wp_o = {(-1)}^{(o/2)+1}$ (85)

The degree damping of damping is modulated by the $\kappa_{do}$ coefficient. When this coefficient is small, noise may develop in the divergence fields. In contrast, when large values of $\kappa_{do}$ are used, the vertical velocity field can be damped too much resulting in such problems as significantly lowered precipitation maximums etc., so a compromise must be made in order to suppress small grid scale noise in the divergence fields while avoiding excessive damping. Note that higher order damping acts more like a filter which damps noisy features but leaves the divergence field less damped overall.

As an example, the second order ($o=2$) divergence damping can be expressed as

$\displaystyle D_{u\,d}^2 =$ $\displaystyle {\frac{\kappa_{d\,2}}{\mu_d}} \,
{\frac{\partial}{\partial x}}\left(\nabla \cdot{\bf V}\mu_d\right)$ (86)
$\displaystyle D_{v\,d}^2 =$ $\displaystyle {\frac{\kappa_{d\,2}}{\mu_d}} \,
{\frac{\partial}{\partial y}}\left(\nabla \cdot{\bf V}\mu_d\right)$ (87)

A balance between minimizing the smoothing and including the damping effects (i.e. prescribing $\kappa_{do}$) is achieved through a combination of using values from the literature and testing. For the example of second order damping (Eq.s 87-88), Shapiro (1991) tested values of $\kappa_{d2}$ in the range 2$x$10$^5$ to 6$x$10$^6$ m$^2$ s$^{-1}$ derived assuming a grid spacing $d=100$ km. For the current model, the value $k_{d}$= 0.028 was used with explicit forward-in time differencing of the prognostic momentum equations, which results in $\kappa_{d2}$= 9.3$x$10$^5$ m$^2$ s$^{-1}$ for $d=100$ km (assuming a Courant number of 1), which obviously falls about in the middle of values that Shapiro tested. Model results can be quite sensitive to divergence damping when an explicit time integration method for the momentum equations is used, depending on the chosen value for $k_{d}$. Slightly larger-than-optimum values damp the vertical velocity fields quite a bit in time, and too low values can result in numerical instabilities. Even reasonable values will tend to dampen the vertical velocity fields such that the amplitudes decrease noticeably in time, but again, use of higher order derivatives can greatly alleviate this problem at fairly low cost.

A forward-backward time-splitting scheme for treating gravity waves is used to integrate the horizontal momentum equations which reduces the need for divergence damping, so that a lower value for $k_{d}$ can be used. For example, Skamarock et al. (2005) use a value of $k_{d}$=0.01 for the WRF model using second order damping. Testing using this value together with the forward-backward scheme was found to give good results, however there was still a general damping over time of wavelengths larger than $2d$. For this reason, a very high order damping is currently used ($o$=4) in order to more selectively damp small scale noise while leaving the overall divergence field relatively intact. The fourth order divergence damping terms can be expressed as

$\displaystyle D_{u\,d}^4 =$ $\displaystyle - {\frac{\kappa_{d\,4}}{\mu_d}} \,
{\frac{\partial^3 }{\partial x^3}}\left(\nabla \cdot{\bf V}\mu_d\right)$ (88)
$\displaystyle D_{v\,d}^4 =$ $\displaystyle - {\frac{\kappa_{d\,4}}{\mu_d}} \,
{\frac{\partial^3 }{\partial y^3}}\left(\nabla \cdot{\bf V}\mu_d\right)$ (89)

The use of fourth order damping maintains a very dynamic pattern with good structure in the vertical velocity field during the entire integration period. See the section on time discretization for more details on the forward-backward time-split integration method.