Evaluating and Improving Wind Forecasts over South China: The Role of Orographic Parameterization in the GRAPES Model

The Global/Regional Assimilation and Prediction System for the Tropical Mesoscale Model (GRAPES_TMM) has presented a high wind speed bias over South China (Zhong and Chen, 2015). The wind bias causes deficiencies in model performance——for example, an underestimation of precipitation and an overestimation of temperature over South China (Zhong et al., 2016), and an overestimation of snowfall over East Asia (Choi and Hong, 2015). Although the bias was alleviated by the inclusion of orographic drag parameterization (Zhong et al., 2016; Choi et al., 2017), it persists in the lower troposphere, in particular over complex terrain and in the planetary boundary layer (PBL). High wind bias is a widespread phenomenon over mountains and valleys (Cheng and Steenburgh, 2005). Nevertheless, wind speed bias also exists in other mesoscale models [e.g., the Weather Research and Forecasting (WRF) model (Skamarock et al., 2008; Lorente-Plazas et al., 2015)].

The orographic gravity wave drag (GWD) can exert forces by gravity wave dissipation, whereas the body forces on the flow can either accelerate or decelerate atmospheric winds (Kim et al., 2003). The wind bias thus may be caused by an unrealistic simulation of the drag generated by the subgrid-scale gravity wave (Lindzen, 1981; Matsuno, 1982; Miller et al., 1989). Thus, GWD parameterization is generally recognized as the critical component for most models (Kim and Arakawa, 1995; McLandress et al., 2012; Choi and Hong, 2015). The topography adopted in GRAPES_TMM is derived from horizontal interpolation by a four-point averaged scheme, which may cause the smoother topography used in the model. It probably leads to an underestimation of the total orographic drag and causes high wind speed bias, especially in the lower troposphere over South China. It is argued that unresolved topographic effects (UTEs) produce an additional drag to that generated by vegetation, which leads to an overestimation of the wind speed in WRF [Jiménez and Dudhia, 2012(JD12)]. The influences of UTEs include small-scale orography (SSO) effects (e.g., SSO drags derived by several hills inside a horizontal grid). This subgrid orographic drag was named turbulent orographic drag (Belcher and Wood, 1996), and it may be of the same order of magnitude as the GWD (Sandu et al., 2016). (Wood et al., 2001) represented the turbulent orographic drag using an explicit orographic stress profile, which considers the effects of stability instead of the effective roughness length (ERL). (Belcher and Wood, 1996) argued that a low-order closure model produces larger orographic drag because it overestimates the shear stress, where drag is related to turbulence.

Several schemes of unresolved orography have been developed to alleviate wind speed bias (Georgelin et al., 2000; Beljaars et al., 2004; Rontu, 2006; Sandu et al., 2016). Among them, three basic concepts are adopted in the parameterization of UTEs. The first is to use an ERL (Fiedler and Panofsky, 1972), which enhances the grid cell vegetative roughness length; the second proposes the introduction of the sink term in momentum equations (Wood et al., 2001; Wilson, 2002); and the third is based on the momentum sink term (MST) method, which multiplies friction velocity according to the terrain characteristics (Jiménez and Dudhia, 2012). These three schemes can improve the performance of numerical weather prediction (NWP) models, especially on surface wind simulations (Milton and Wilson, 1996; Rontu, 2006; Jiménez and Dudhia, 2013; Gòmez-Navarro et al., 2015; Lee et al., 2015; Sandu et al., 2016). Parameterizations based on the MST method show certain advantages over those adopting the ERL method. The MST method shows no limitation of the lowest model level and takes orographic features into account (Beljaars et al., 2004; Jiménez and Dudhia, 2012). However, the smoother topography used in the model can potentially lead to an underestimation of the SSO drags. For instance, the height of the mountain used in the model is lower than reality, especially over small hills within a horizontal gird. It has been argued that wind speed overestimation at mountainous regions would be due to lack of representation of SSO drag in NWP models [e.g., in the WRF model (Jiménez and Dudhia, 2012)].

Another related characterization parameter of the surface wind in NWP models is friction velocity. It was found that excessive convection-induced turbulent mixing under free convection can lead to an overestimation of friction velocity and impose excessive surface friction (Liu et al., 2004). The overestimation of friction velocity can thus reduce surface momentum and lead to excessively weak surface winds. Underprediction of daytime surface winds results from an underestimation of downward momentum flux (Zhang and Zheng, 2004) or overestimation of surface stress (Liu et al., 2004). (Lorente-Plazas et al., 2016) enhanced the friction velocity by considering the turbulent orographic drag effect and stability effects, and the results demonstrated that the inclusion of stability effects improved simulated surface winds by alleviating the systematic daytime underestimation of the original scheme.

The present work aims to alleviate wind speed bias in the GRAPES_TMM model. The SSO drag effects on surface winds are parameterized by considering the orographic features based on the JD12 scheme. The SSO drag effects are represented by adding a sink term in the momentum equations with the maximum height of the hills or mountains within the grid box. The modification could be important by compensating for the drag that is underestimated by the smoother terrain used in the static data. The modified scheme is implemented and coupled with the PBL parameterization scheme within the model physics package. The results are verified against wind observations over South China and reanalysis data.

The model used in this study is the new-generation Mesoscale Atmospheric Regional Model System (MARS) in operation over South China based on GRAPES_TMM (Zhong et al., 2016). MARS focuses on short-range forecasting of mesoscale convective systems over South China. The resolution was updated to within 3 km in 2016, and it provides data four times daily (initialized at 0000, 0600, 1200 and 1800 UTC). The physics package used in this paper is the same as that used in (Zhong and Chen, 2015), but without the convective cumulus scheme. The simulation domain comprises 913× 513 grid points, with a horizontal resolution of 3 km, as shown in Fig. 1a. There are 60 layers in the vertical direction, and both the initial and lateral boundary fields are obtained from the 0.125^°× 0.125^° forecast fields of the European Centre for Medium-Range Weather Forecasts (ECMWF), with the lateral boundary fields updated every 6 h.

MARS uses a modified orographic parameterization scheme by including mountain blocking drag based on (Kim and Arakawa, 1995). This study focuses on the alleviation of wind bias through the implementation of the SSO parameterization (SSOP) scheme in MARS, based on the momentum-conservation equation (Jiménez and Dudhia, 2012). The topographic height (h) used in MARS is derived from horizontal interpolation by a four-point averaged scheme. In this paper, correction of the non-dimensional parameter of the topography ∆²h_i,j is updated by using the maximum orography height (h') from the U.S. Geological Survey at a horizontal resolution of 30 arc-seconds (Gesch et al., 2002). The differences between h' and h are shown in Fig. 1b. It can be seen that the differences in subgrid orography can reach more than 200 m around Xinyi, west of Guangdong Province, South China. Thus, it employs the maximum height h' of the mountains as compensation for drag being underestimated by the smoother terrain h used in the static data: \begin{equation} \Delta^2h'_{i,j}=0.25(h'_{i+1,j}+h'_{i-1,j}+h'_{i,j+1}+h'_{i,j-1}-4h'_{i,j}) , \ \ (1)\end{equation} where positive values of ∆²h represent valleys and negative values indicate mountains; near-zero values indicate plains. The effects of the unresolved topography are parameterized, introducing a corresponding parameter c_t as a modulation of the surface drag associated with topography height in the calculation of the tridiagonal matrix elements (Hong et al., 2006) for momentum tendencies in PBL parameterization: \begin{equation} {AM}_{ v}={F}_{ sd} .\ \ (2) \end{equation}

Figure 1.  (a) Domain of MARS (3 km). Color shading denotes the orography, and the red frame shows the area of concern in (b). (b) The differences between h' and h; the dots denote the observation sites.

The vector of momentum tendency M _v is solved using the matrix A and forcing F _sd, where the main diagonal of A on the first model level is defined as A₁: \begin{eqnarray} A_1&=&(1+\delta_1) ;\ \ (3)\\ \delta_1&=&2c_{ t}\rho g\dfrac{u_*^2}{s\Delta z}dt .\ \ (4) \end{eqnarray} Here, the subscript "1" stands for the first model level, ρ is the air density, g is the acceleration due to gravity, δ₁ represents the sink term, S is the wind speed at the first model level, u_* is the friction velocity, and ∆ z is the thickness of the first model layer. c _t is a function of ∆²h'_i,j and the standard deviation σ' _sso: \begin{equation} c_{ t}=\left\{ \begin{array}{l} 1,\;\Delta^2h'_{i,j}>-20,\;\sigma'_{ sso}<e\\ \ln\sigma'_{ sso},\;-10<\Delta^2h'_{i,j},\;\sigma'_{ sso}>e\\ \alpha\ln\sigma'_{ sso}+(1-\alpha'),\;-10>\Delta^2h'_{i,j}>-20,\;\sigma'_{ sso}>e\\[1mm] \dfrac{\Delta^2h'+30}{10},\;-20>\Delta^2h'_{i,j}>-30\\[3mm] 0,\;\Delta^2h'_{i,j}<-30 \end{array} \right. ,\ \ (5) \end{equation} where α'=∆²h'+20/20 and e is the natural logarithm base. The standard deviation σ' _sso is also derived from h': \begin{eqnarray} \sigma'_{ sso}=\sqrt{\sum\dfrac{(h'-\overline{h'})^2}{N-1}} , \ \ (6)\end{eqnarray} where \(\overline{h'}\) is the mean value in the grid box and N is the number of the grid box. A comprehensive description of the scheme and its numerical discretization is presented in the Appendix.

The scheme is tested against wind observations over South China——a complex terrain region that includes the Nanling Mountains and Yunwu and Yunkai Mountains (YYMs), as well as the longest coastline in China (Fig. 2a) and a complex underlying surface (mountains, plains, valleys, lakes, and so on). May 2016 is selected as the verification period; the main outbreak of the South Asian monsoon occurs during this period. The strong southerly winds bring abundant moisture and cause severe torrential rain, especially given the large-scale and small-scale mountains. Xinyi is characterized by a trumpet-shaped topography over the west of Guangdong Province (Fig. 2b), where the north is dominated by the northern and western steep ridges of the YYMs and the south is covered by low valleys and areas of flat terrain.

A total of 42 observational stations with conventional meteorological element measurements are used in the zoomed area in Fig. 2, and around 8000 stations in the whole domain of MARS (Fig. 1). These sites over the YYMs cover most of the region with a fine network. For instance, there are fifteen stations located over the flat region and five stations over the mountain region, as well as fourteen stations in the valleys and eight in the flat areas between mountains and valleys. The wind coverage shows that it is dominated by southeasterly wind over the flat areas in the south of Xinyi, with wind speeds of around 4 m s^-1. Torrential rain occurred in association with the circulation over Xinyi on 20 May 2016. The 3D coverage of the topography, 24-h rainfall and wind distribution are shown in Fig. 3. It can be seen that precipitation was concentrated over the trumpet-shaped topography of Xinyi. The 24-h accumulated rainfall reached 464 mm, with southeasterly winds over the south and southeast of the Yunwu Mountains.

Figure 2.  (a) The topography around Xinyi in the west of Guangdong Province and (b) the ∆²h' (negative values represent mountains, and positive values are valleys). The wind vectors (units: m s^-1) in (a) are the five-day wind observation averages from 0800-1200 UTC, 16-18 May 2016. The red contours in (a) denote surface wind speed (units: m s^-1).

Figure 3.  The 3D topographic features of Xinyi (color-shaded; units: m), 24-h rainfall at 1200 UTC 20 May (contours; units: mm), and surface wind distribution at 0000 UTC 20 May 2016 (vectors; units: m s^-1).

The main diagonal A at the first model level in Eq. (3) is related to the wind speed, friction velocity, thickness of the first model layer and topographic features c _t. For a given forcing term in Eq. (2), the wind speed becomes smaller as A₁ gets larger, and vice versa. Also, the value of A₁ is proportional to the magnitude of c _t and friction velocity. Figure 4 shows the distribution of c _t and A₁ for the 50-time-step integration. It can be seen that the distribution of c _t and A₁ shows good consistency with the topographic features in Fig. 1. The distribution of c _t in the steep mountainous region, especially over the Tibetan Plateau, is larger than those over the flat regions. Thus, the scheme could alleviate the wind speed in these mountainous regions as momentum tendency becomes smaller after considering subgrid orographic parameterization.

Figure 4.  The distribution of (a) topographic function c _t and (b) sink term A₁ for the 50 time-step integration.

Figure 5.  Comparison of the friction velocity (shaded; units: m s^-1) of the 24-h simulation between (a) CTL and (b) ORO.

In addition, friction velocity is also an impact factor with respect to A₁ according to Eq. (4). Three experiments are set up in this study. Two of these experiments are the ORO and CTL experiments, which represent the SSOP scheme and control simulation, respectively. Figure 5 shows the feedback of the SSOP to the friction velocity in these two experiments. It can be seen that the friction velocity is alleviated after using the SSOP scheme in the ORO experiment. As enhanced friction velocity can lead to a larger A₁ in Eq. (4), and thus weaken the momentum tendency, this could lead to an underestimation of daytime winds (Liu et al., 2004). In this study, the results show that the friction velocity also becomes smaller, since the overestimated wind speed is alleviated (see section 3.2). However, the JD12 scheme may cause the underestimation of the daytime wind speed. Some other measures should be considered to alleviate this phenomenon; for instance, the effects of atmospheric stability could be included in the SSOP scheme (Lorente-Plazas et al., 2016), which is not discussed in this investigation.

This study focuses on the subgrid orographic drag effects on low-level winds, especially the landing southerly winds from the South China Sea (e.g., southwesterly winds brought by the monsoon), which have a strong influence on the formation of the torrential rain over this complex terrain (Lin, 2006). As mentioned, this investigation includes three experiments, including the ORO and CTL experiments introduced above, but also a third experiment of the JD12 scheme with the ORO experiment to examine the subgrid orographic effects. The Medium-Range Forecast (Hong et al., 2006) PBL scheme is used to verify the effects from the SSOP scheme, which is an independent code in the model physics process.

The surface winds of the ORO and CTL experiments are compared with observations in Fig. 6. The two experiments share similar simulated circulations, whereas the CTL experiment exhibits an obvious overestimation of wind speed. The magnitude of the wind speed over the flat regions reaches more than 3 m s^-1 in the CTL experiment, while the observations are around 1 m s^-1. More overestimations are found over the mountains. For instance, the simulated winds reach more than 5 m s^-1 over the east and northwest of Xinyi, with observations around 3 m s^-1. In contrast, the ORO experiment provides a much better simulated result than the CTL experiment, through a consistent magnitude with the observation over the flat regions and a slightly overestimated wind speed over the mountains.

The 850-hPa wind simulations are compared between the two experiments in Fig. 7. Although both experiments share a similar simulation of the low-pressure vortex in the west of Xinyi, the southerly winds at the eastern edge of the vortex exhibit huge differences. The CTL experiment shows higher wind speeds than in the ORO experiment, with a smaller southerly wind zone over the eastern edge of the vortex. Another apparent difference in the two experiments is the southerly wind of the coastal regions over the southwest of Guangdong and the Beibu Gulf. The ORO experiment shows a larger southerly wind over the coastal areas, whereas the abundance of moisture accompanying the southerly wind may have a substantial influence on the torrential rain over the terrestrial mountainous regions in the north.

Figure 6.  Comparison of the surface wind between the 12-h simulation (color-shaded) and the observation [colored dots (same color scale bar as the forecast); units: m s^-1] at 0000 UTC 20 May, in which the contours represent the topography, for (a) CTL and (b) ORO.

The verification of the 850-hPa wind simulations over the zoomed area in Fig. 2a by the two experiments is analyzed against reanalysis data in Fig. 8. Both experiments underestimate the winds over the Yunkai Mountains and northwest mountainous regions. The CTL experiment overestimates the wind over the south and east of Xinyi, which means the simulated southerly wind is too strong, especially over the flat regions over the south of Xinyi and over the Yunwu Mountains in the east of Xinyi. The ORO experiment provides much better simulated southerly winds over these regions. The simulated winds show a slight overestimation in some parts of the flat regions and an underestimation over the Yunwu Mountains.

Figure 7.  Comparison of the 850-hPa wind speed (color-shaded; units: m s^-1) and relative vorticity (contours; >0.4×10^-3 s^-1) for the 24-h simulation at 1200 UTC 20 May between (a) CTL, (b) ORO and (c) ORO minus CTL, and the black frame shows the location of the vortex.

Figure 8.  As in Fig. 7 but a comparison of the 850-hPa wind speed against the analysis (units: m s^-1): (a) analysis minus CTL; (b) analysis minus ORO.

Figure 9.  Comparison of the 24-h simulated wind profile averaged (a) over the whole domain and (b) in the zoomed area in Fig. 2.

Figure 9a shows the vertical profile of domain-averaged 24-h simulated wind speed. The two schemes have no significant influence on the wind speed simulation, especially at layers higher than 800 hPa. Both schemes alleviate the wind speed between 1000 hPa and 800 hPa, and the 24-h wind speed simulation by the ORO scheme is slightly smaller than that of the JD12 scheme. The average wind profiles of the two experiments and the JD12 scheme over the coastal regions of southwestern Guangdong are compared in Fig. 9b. The influence of the scheme is more complicated, since it is related to the co-effects of the mountainous areas and coastal regions, and the feedbacks of the scheme to the wind speed can nearly reach 300 hPa, which is probably caused by the vertical transport of momentum and vertical diffusion in the PBL. Though both schemes can alleviate the wind speed, the ORO scheme provides a weaker wind speed than the JD12 scheme. Note, however, that revealing the sensitivity of the differences to the chosen scheme requires more experimentation, which is beyond the scope of this study.

A monthly verification is executed to examine the stability and monthly performance of the SSOP scheme. This experiment is initialized twice daily at 0000 and 1200 UTC in May 2016, and the results are verified by the analysis field and observational station data, including the wind speed at 1000 hPa, 925 hPa and 850 hPa by the analysis field, and the 10-m wind (V10m) and 2-m temperature (T2m) by the observations. Parameter verification is conducted by calculating the root-mean-square error (RMSE) between the forecast departure and real-case departure, the calculation of which is as follows: \begin{equation} { RMSE}=\left[\dfrac{1}{N}\sum(S-A_{ v})^2\right]^{\frac{1}{2}} , \ \ (7)\end{equation} where S is the forecast value, A _v is the analysis value (observation data), and N is the number of grid points in the verification region. The verification of wind speed at the isobaric level is employed in the ECMWF analysis field as the analysis value. The observations used in the verification of V10m and T2m are the observed ground-level surface synoptic observations (SYNOP) data collected by the South China Regional Meteorological Center at about 10 000 observation stations, some of which are shown in Fig. 1b.

As can be seen from Table 1, the average 24-h and 48-h wind speed bias is greatly alleviated in the ORO experiment. For instance, the wind forecast error at 1000 hPa is alleviated from 4.67 m s^-1 to 2.18 m s^-1, which is also verified using the V10m observations, with an error reduction of 2.99 m s^-1 and 3.1 m s^-1 from the 24-h and 48-h forecast (Table 2). The reduction in 24-h wind speed RMSE at 925 and 850 hPa is 1.05 and 0.34 m s^-1, respectively; plus, there are 48-h alleviations of 0.79 m s^-1 and 0.15 m s^-1. The magnitude of the wind speed error reduction becomes weaker with a longer integration. On the other hand, the error reduction in the lower layer is more significant than that in the higher level, with a slight improvement in the T2m forecast (Table 2). The overestimated wind speed is thus reduced by including the SSO drags.

Figure 10 compares the monthly averaged surface wind between the 24-h simulation and the observations. It can be seen that the simulation without the SSOP scheme overestimates the surface wind over most of the region, especially over the flat regions in the south of Xinyi. The simulated winds over these regions reach more than 6 m s^-1, while the observations are about 1-2 m s^-1. Although both simulations can capture the strong wind speeds over the mountains, the simulation with the SSOP scheme provides a much better wind distribution over these regions, especially the flat regions over the south of Xinyi. Nevertheless, the simulated wind shows slight overestimation in some parts of the flat regions and underestimation over the Yunwu Mountains.

Figure 10.  As in Fig. 6 but a comparison of the monthly averaged surface wind between the 24-h simulation and the observations (units: m s^-1): (a) CTL; (b) ORO.

Unresolved SSO drags are parameterized to improve the forecasting of surface winds in GRAPES_TMM. The SSO drags are represented by adding a sink term in the momentum equations. A modification of the terrain features ∆²h is calculated by the maximum height of the mountains within the grid box as a compensation for the drag. The scheme is implemented and coupled with the PBL parameterization scheme within the model physics package. A monthly simulation using the modified scheme outperforms the surface wind estimations and the default simulation over the complex terrain areas located in the southwest of Guangdong.

It is found that surface wind speed bias is greatly alleviated by adopting the SSOP scheme; this is also the case for surface temperature bias and wind bias in the lower troposphere, especially in the PBL. The monthly verification results show the wind forecast error at 1000 hPa to be alleviated from 4.67 m s^-1 to 2.18 m s^-1, with an error reduction of 2.99 m s^-1 and 3.1 m s^-1 from the 24-h and 48-h forecast, respectively. The reduction in 24-h wind speed RMSE at 925 and 850 hPa is 1.05 and 0.34 m s^-1, respectively; plus, 48-h alleviations of 0.79 m s^-1 and 0.15 m s^-1 are also found. The verifications between the simulation and observations show that, although both simulations can capture the strong wind speeds over the mountains, the simulation with the SSOP scheme provides a much better wind distribution over these regions, especially over the flat regions in the south of Xinyi. Nevertheless, the SSOP scheme still shows a slight overestimation over some parts of the flat regions and an underestimation over the Yunwu Mountains.

The results presented here show that the SSOP scheme in this paper to some extent outperforms the JD12 scheme. The SSO drags are strengthened by using a more realistic terrain height of the mountains, whereas the drags are underestimated with the current four-point averaging interpolation scheme in the static data. The SSOP scheme in this paper should be independent to the horizontal resolution used. Thus, more verifications should be designed with a more realistic terrain to evaluate SSO drags with different horizontal resolutions, and a more detailed verification of the performances with respect to the diurnal wind speed over other complex terrains should be examined. Moreover, the need to include stability effects should be considered to improve the diurnal characteristics of wind speed representation during the daytime and at nighttime.

Acknowledgements. Special thanks are given to the editors for the formula normalization. We also thank the reviewers for their helpful comments. This study was supported by the National Natural Science Foundation of China (Grant Nos. 41505084, 41275053 and 41461164006), the China Meteorological Administration Special Public Welfare Research Fund (Grant Nos. GYHY201406003 and GYHY201406009), the Guangdong Meteorological Service Project (Grant No. 2015B01), and the Guangdong Province Public Welfare Research and Capacity Construction Project (Grant No. 2017B020218003).