Filtering

Finally, in order to find an optimum level of CPU times which permits the use of a grid spacing which is larger than the minimum value on the grid near the Poles but avoids numerical instabilities arising from the decreasing grid spacing, Fourier Filtering is used above some critical latitude, $\phi_{r,crit}$. We use the classic filter weight, $w_f$, form (Haltiner and Williams (1980)) applied to the Fourier coefficients:

$$w_f(k) = {\left[ {\frac{{\rm cos}(\phi_r)}{{\rm cos}(\phi_{r,crit... ...{1}{{\rm sin}(\pi k/n)}} \right]}^m \hskip1.in \left(0\leq w_f(k) \leq 1\right)$$ (314)

where $k$ is the wave number and $n$ is the number of grid points along the latitude $\phi_r$. Finally, $m$ is a parameter which controls the filter strength. It can range from 1 (minimum) to 2 for increased stability (smoothing), however, the default value herein is to use $m=1$ in order to minimize the smoothing while ensuring numerical stability.

The variables which are filtered are $p_s$ and $\mu\phi$, where $\phi$ represents $u$, $v$, and $T$. Filtering is performed for $p_s$, $\mu \, u$, $\mu \, v$ and $\mu \, T$ after each time split within the RK loop.

It has been seen that for coarse resolution tests at very high latitudes where steep topographic gradients occur, an improved filtering of $p_s$ can be done by filtering the perturbation of surface pressure defined as

$$p_s^\prime = p_s - p_{s\,std}$$ (315)

where the time constant standard surface pressure is defined using the WRFv3 definition as

$$p_{s\,std} = p_0 \, {\rm exp} \Bigg\lbrace - {\frac{T_0}{A_T}} \,... ...A_T}} \right)}^{2} \,-\, {\frac{2 \Phi_s}{A_T R_d}} \right]}^{1/2} \Bigg\rbrace$$ (316)

where $A_T=50$ K, $p_0=1000$ HPa and $T_0=288$ K. Once $p_s^\prime$ has been filtered, $p_{s\,std}$ is added back to get an updated $p_s$. Thus, mass is conserved and the surface pressure fields are more effectively smoothed (again, for the special conditions mentioned above). The same method is used for the potential temperature using a background time constant temperature, which is defined as

$$T_{s\,std} = T_0 \,+\, A_T {\rm ln}\left( {\frac{A_\eta + B_\eta p_{s\,std} }{p_{s\,std}}} \right)$$ (317)

Thus,the filtering can be expressed as

$${\overline{\pi}_s} = {\cal F}\left( \pi_s - p_{s\,std} \right) \,+\, p_{s\,std}$$ (318)

$${\overline{\theta_\rho}} = {\frac {\overline{\theta_\rho\mu}}{\overline{\mu}}} = {\frac{1}{{... ...eft( \theta_\rho\mu - T_{s\,std}\mu/\Pi \right) \,+\, T_{s\,std}\mu/\Pi \right]$$ (319)

$${\overline{u}} = {\frac {\overline{u\mu}}{\overline{\mu}}} = {\frac{1}{{\overline{\mu}}}} {\cal F}\left( u\mu\right)$$ (320)

$${\overline{v}} = {\frac {\overline{v\mu}}{\overline{\mu}}} = {\frac{1}{{\overline{\mu}}}} {\cal F}\left( v\mu\right)$$ (321)

$${\overline{w}} = {\frac {\overline{w\mu}}{\overline{\mu}}} = {\frac{1}{{\overline{\mu}}}} {\cal F}\left( w\mu\right)$$ (322)

$${\overline{\alpha}} = {\frac {\overline{\alpha\mu}}{\overline{\mu}}} = {\frac{1}{{\over... ...F}\left( {\mu\alpha}- {\mu\alpha}_{std} \right) \,+\, {\mu\alpha}_{std} \right]$$ (323)

$${\overline{\Phi}} = {\cal F}\left( {\Phi}- {\Phi}_{std} \right) \,+\, {\Phi}_{std}$$ (324)

where the ${\overline{\phi}}$ is used to represent the filtered variable, $\phi$, and ${\cal F}$ represents the Fourier Filter operator. ${\overline{\mu}}$ represents $\mu$ computed from the filtered surface pressure, ${\overline{\pi_s}}$. These mass variables are conserved during the filtering. $\alpha$ is computed from Eq. 258 using ${\overline w}$, then the non-hydrostatic pressure, $p$, is computed using Eq. 259 using ${\overline{\alpha}}$ (these two variables are not filtered if the non-hydrostatic option is off). Note that for higher resolutions, the need for the background state for $\theta$ and $\alpha$ in the filtering diminishes. Finally, either Eq. 325 or Eq. 326 is used, depending on whether $\Phi$ or $\mu\alpha$ is used as a prognostic variable within the non-hydrostatic module.

From Eq. 322, it can be seen that $u\mu$ (and equivalently $U$) is filtered, thus the horizontal mass flux, $F$, has been filtered. For the scalar water substance equations and the non-hydrostatic continuity equation, $F$ only appears in the longitudinal advection terms, and numerical stability can be maintained simply be using the filtered $F$ value. Therefore, the scalar water substance variables and geopotential are not filtered, thus the scheme is quite economical.

Finally, note that since the RK time step is smaller than $\Delta t$ for all RK steps before the final one, the critical latitude for maintaining numerical stability is a maximum for the first RK step, then decreases for each RK step (as to maintain numerical stability). This saves considerable CPUs as the Fourier transforms are computed over a larger range of latitudes as the RK stepping progresses. For example, for $\phi_{r,crit}$ 45$^o$, we compute the model time step, $\Delta t$, based on the grid distance at this latitude. This corresponds to the value used in the last RK step. For the case using using RK3, the critical latitudes for RK(step 1) = 74$^o$ and RK (step 2) = 66$^o$ (these values are slightly lower than the maximum theoretical values in order to provide a safety margin). Note that filtering is only applied during dynamics computations (it is not applied after micro-physics and turbulence computations, which are performed after the dynamics). An example forecast on the global grid is shown in Fig. 10.

Figure 10: A sample forecast from ASP on the global grid.