A Review on GRAPES-TMM Operational Model System at Guangzhou Regional Meteorological Center

The Global/Regional Assimilation and Prediction System-Tropical Monsoon Model (GRAPES-TMM ^[1]) has been in operation since 2005 at the Guangzhou Regional Meteorological Center (GRMC). The earliest version of the model systems based on GRAPES-TMM included the Guangzhou large-scale model for tropical and marine meteorology (with a horizontal resolution of 36 km and 31 vertical layers) and Guangzhou mesoscale model (12 km, 31 layers). In 2011, the first version of the Tropical Regional Atmospheric Model System for the South China Sea (TRAMS-v1.0, 36 km) and the Mesoscale Atmospheric Regional Model System (MARS, 12 km) were successfully developed based on the GRAPES-TMM, focusing on the theoretical research and operational application of numerical weather prediction (NWP) in the tropical region.

With the increasing computing capabilities of computers, the NWP systems in the GRMC have experienced a period of rapid development since the establishment of the Key Laboratory of Regional Numerical Weather Prediction (LARP) through the province-ministry cooperation between Guangdong Province and China Meteorological Administration in December 2012. After more than one year of experimental evaluation in 2015, the TRAMS-v1.0 was upgraded to v2.0 in June 2016. TRAMS-v2.0 adopts a higher horizontal resolution of 18 km and 65 vertical layers (Chen et al. ^[1]). The related improvement of model dynamics and physics includes the development of a three-dimensional (3D) reference atmosphere technique, an improvement on the calculation accuracy of nonlinear terms, coupling technique for model dynamics and physics (Chen et al. ^[2]), the parameterization of gravity wave draw induced by sub-grid orography (GWDO, Zhong et al. ^[3]), the improvement of cumulus parameterization (Xu et al. ^[4]) and the coupling technique for convective process and microphysical process (Xu et al. ^[5]). Evaluations indicate that TRAMS-v2.0 shows an overall improved performance on typhoon forecasting than the v1.0 did.

In the meantime, the rapid update cycling (RUC) forecasting system (1km) has been in operation based on the "Tianhe-2" Super Computer in Guangzhou since 2015. The RUC-1km, focusing on the short-range weather forecasting mainly in Guangdong Province, adopts several model techniques as those in TRAMS-v2.0. The w_damping technique is adopted to increase the stability of the model, improve the water vapor transport scheme and integrate the cloud analysis system. In December 2018, MARS-v1.0 was subsequently upgraded to 2.0, which adopted a horizontal resolution of 3 km (Zhong et al. ^[6]). MARS-v2.0 adopts similar techniques as those in TRAMS-v2.0 and RUC-1km. Besides, it integrates the improved parallel nudging technology, the radar-upgraded cloud analysis system, and the enhanced calculation accuracy of nonlinear terms.

In 2019, both TRAMS and MARS expanded the forecasting domain in response to operational needs and the "Belt and Road Initiative" and the "21^st Century Maritime Silk Road Initiative". In July 2019, TRAMS-v2.0 was upgraded to v3.0, which adopts a horizontal resolution of 9 km. TRAMS-v3.0 had an optimized input/output (I/O) module and model initialization, which effectively improves the parallel efficiency of the model. Furthermore, the dynamic framework of the model has also been optimized. Fig. 1 shows the current forecasting system in the GRMC, including the TRAMS-v3.0 (9km), MARS-v2.0 (3km), and RUC (1km), which focuses on the typhoon forecasting in the west of the Pacific (South China Sea), operational NWP in southern China and short-range weather forecasting in Guangdong Province, respectively.

Figure 1.  The domain of TRAMS-9km, MARS-3km and RUC-1km (9-3-1).

This review summarizes the related progress in the field of numerical atmospheric modeling based on GRAPES-TMM in the past two decades. The paper is organized as follows. Section 2 and 3 synthesize advances in the field of model dynamical core and physical process. Section 4 reviews the progress in model evaluation. The conclusions are given in section 5.

With the increase of global model resolution and improvement of the numerical model technique, the focus of the regional numerical weather model is turning to convection-scale weather forecasting, which requires the model to be a non-hydrostatic model. A high-resolution regional model it includes Rossby wave, gravity wave and acoustic wave (Kasahara and Qian ^[8]). Rossby wave is a slow synoptic-scale perturbation, and the vertical perturbation of gravity wave and the acoustic wave is significant (Daley ^[9]). Nevertheless, the propagation speed of external gravity wave and acoustic wave is very fast, which has substantial impacts on computational precision and thus may lead to computational instability. Therefore, it requires a shorter time step and a higher precision discretization scheme.

Reference atmosphere is one of the useful methods to linearize the model equations. The dynamical core of the previous version of GRAPES-TMM used a longitude-latitude horizontal grid with a 1D profile of the hydrostatic atmosphere. The new version adopts a 3D reference atmosphere scheme in which the reference atmosphere is not only a function of height but also a function of longitude and latitude. With the introduction of the 3D reference scheme, the atmosphere is divided into the basic state which satisfies the following hydrostatic equation and perturbation which is non-hydrostatic:

where λ(ϕ) denotes the longitude (latitude) of spherical coordinates and z denotes the vertical coordinates. $\bar{\Pi}, \bar{\theta} $ represents the basic state of reference atmosphere which satisfies the hydrostatic equation. Π' θ' represents the perturbation deviating from the basic state. Unlike the large magnitude of the forecasted perturbation in the model derived from the 1D reference, the magnitude of the forecasting perturbations is significantly reduced with the adoption of the 3D reference method (Fig. 2).

Figure 2.  Comparison of initial potential temperature perturbation for 1D (red) and 3D (blue) reference atmosphere (x-axis represents potential temperature perturbation; y-axis represents model level).

In May 2015, a one-month experiment was performed to examine the difference of mean absolute errors between 1D and 3D schemes. Fig. 2 shows the mean absolute error of geo-potential height and temperature at different levels. It can be seen that the 3D scheme exhibited a much smaller error of temperature and geo-potential height than the 1D scheme do. The mean absolute error of temperature and geo-potential height for the 1D scheme ranged from 0.8 to 2.1K and 5.6 to 15.6 gpm at different levels. Nevertheless, the 3D scheme showed a general smaller error ranging from 0.4 to 1.2K and 2.7 to 7.3 gpm. It should be noted that the 3D reference atmosphere scheme was first developed and implemented in the GRAPES model at the Guangzhou Regional Meteorological Center by Chen et al. ^[10]. Due to its extensive applicability and reliable performances, the scheme was thereafter implemented into the GRAPES_GFS (GRAPES, Global Forecast System) in the National Meteorological Center in Beijing (Su et al. ^[11])

The dynamics and physical processes are two key technologies in the development of NWP models. The dynamics are calculated on the model grids, including advection, adjustments, and diffusions. The physical processes are mainly related to the phase-transition on sub-grid processes which cannot be solved directly, including planetary boundary layer (PBL) process, tropical cumulus convection, topographic parameterization (Zhong et al. ^[6]), and radiations. The model physics is generally derived by the dynamics, and the evolution of dynamics is also affected by physics. Therefore, the effective coupling between dynamics and model physics are inevitable.

To explore the effects of the coupling between dynamics and physics on model accuracy, Chen et al. ^[2] developed a coupling scheme, which includes the feedbacks of temperature and water vapor to Π. Furthermore, the physical feedbacks are calculated in the Helmholtz equation of the model as the implicit solution. The experiment results showed that the coupling scheme improved the overall performances of the model (Fig. 4), especially for the prediction accuracy of the moving speed of typhoons.

Figure 3.  Comparisons of mean absolute error (AE) of 24-h simulated geo-potential height and temperature at different levels between 1D and 3D schemes.

Figure 4.  Comparisons of typhoon track error between the ctl (no coupling scheme) and coupling experiment.

GRAPES-TMM uses the semi-implicit and semi-Lagrangian (SISL) time-difference scheme to improve the computational stability and accuracy of the model, which needs to solve the complex implicit equations. First, the equations need linearization separation, as the nonlinear term has not been included in the implicit solution, and the value of the nonlinear term is usually estimated before the implicit solution. For two steps of the SISL scheme (n, n + 1), the nonlinear term at step n could be resolved before solving implicit equations, and the nonlinear term at step n + 1 remains unknown which needs to be estimated before solving implicit equations. In the original scheme of GRAPES, the nonlinear term at step n + 1 was approximately calculated by the extrapolation method.

where $ \tilde{N}$ is the estimation at step n + 1 and N^f is the time smoothing value. This algorithm is simplified and convenient, whereas accuracy is compromised, especially when the weather system generates and disappears faster.

To improve the accuracy of the model, TRAMS uses a new scheme of nonlinear terms by fractional steps (Chen et al. ^[1]). This algorithm first gets the pre-solution of variables at step n + 1 before solving the Helmholtz equation, then calculates the nonlinear term and physical feedbacks at step n + 1. Therefore, the calculations of physical processes at n and n + 1 could be finished before solving implicit equations, and the physical feedbacks does not need to be directly added to the dynamic process, but as an increment of Helmholtz's right-hand term, which participates in the calculation of Helmholtz equation. By solving the Helmholtz equation, the value of π (intermediate value of barometric pressure) at step n + 1 is obtained, and then the forecasting values of u, v, w and θ at step n +1 are calculated.

The water vapor advection scheme is essential to weather forecasting especially for the tropical regions (e. g., the strong convection of water vapor in a typhoon). Due to the limitations of computational stability by Courant-Fredrichs-Lewy (CFL), the traditional Euler explicit difference scheme must adopt small time steps, which has greatly limited the development of high-resolution NWP models. Therefore, most of the models used a semi-implicit and semi-Lagrangian scheme, which could use a longer time step and theoretically almost had no limitations to the CFL conditions. However, in real settings, due to the limitations of computational accuracy and parallel computing, the time step is much longer than the traditional Euler explicit difference scheme.

GRAPES-TMM uses the quasi-monotone semi-Lagarangain (QMSL) and a second-order moment conservation scheme (Liao ^[12]). Operational evaluations show that this scheme provides a stable performance, but accomplishes with too weak precipitation. Based on the QMSL scheme, the nonlinear constraints used in the linear constraint semi-Lagrangian (LCSL, Pellerin et al ^[13]) scheme are included in the TRAMS (Chen et al ^[14]). For the variable Q at step n + 1, it could be calculated by the Q at step n:

where k represents the grid point value, x (k) and a (k) denote the position vector and displacement from step n to n + 1 at grid k, respectively. x (k) -a (k) denotes the upstream grid position, which is mostly not on the grid point; the grid point value needs to be obtained through linear interpolation by the four grid points around. Experimental simulations show that the modified scheme can improve the forecasting of precipitation and geo-potential height, especially for typhoon forecasting.

The gravity waves (GWs) induced by sub-grid orography has strong impacts on the momentum fluxes in the middle atmosphere (see Kim and Doyle ^[15]; McLandress et al. ^[16]; Zhong et al. ^[7]). Dissipation by GWs could cause synoptic-scale forces known as gravity wave drag (GWD). Parameterization of GWD could improve the overall performances of NWP models by alleviating the systematic wind bias especially in the troposphere (Kim and Arakawa ^[17]; Zhong et al. ^[3]). Moreover, the mountain blocking drag (MBD) also has strong impacts on the wind at the low troposphere. The model including the MBD parameterization shows an improved low-level flow deflection (Lott and Miller ^[18]; Zhong et al. ^[7]).

Zhong et al. extended the GWD scheme based on Kim and Arakawa ^[17] with the inclusion of MBD effects in the GRAPES model. Readers are advised to refer to Kim and Arakawa ^[17] and Zhong et al. ^[3] for a detailed description of the scheme. Here, we present only that part relevant to the GWD scheme by introducing a critical effective height H_c:

where KE represents the kinetic energy, and PE denotes the potential energy. When KE > PE, the atmospheric flow goes over the mountain and triggers orographic gravity waves; if KE < PE, the flow goes around the mountain and generates mountain blocking drag:

where τ₀ is the GWD stress at the reference level and τ_blk is the MBD stress. For each variable in (3) and (4), readers may refer to the introductory part of the orographic drag parameterization by Zhong et al. ^[3]. Evaluations show that the scheme can improve the overall performances on typhoon track forecasting (Fig. 5), especially on the forecasting of typhoon landing.

Figure 5.  Comparisons of typhoon track error between the ctl (no GWDO scheme) and GWDO experiment ^[3].

The cumulus parameterization scheme used in an early version of TRAMS was based on the simplified Arakawa-Schubert Scheme (SAS) by Pan and Wu. ^[19]. This scheme used the quasi-equilibrium assumption and the mass flux concept to adjust the atmospheric temperature, as well as the inclusion of a downdraft scheme that is analogous to the updraft scheme. The scheme was easy to implement in which the entraining and detraining properties of the clouds could be changed to allow experimentation. However, several factors that were essential to convections were neglected, such as the momentum transport and relative humidity feedbacks to the inclusion rate of the cloud-side boundary.

To improve the performances of TRAMS, the New SAS (NSAS) scheme was implemented based on Han et al. ^[20]. The modifications of NSAS includes the following steps:

(1) The NSAS scheme removes unrealistic moisture accumulation in the layer below the inversion, employing the additional diffusion of heat and moisture in the shallow convection (SC) scheme.

(2) The deep convection scheme in SAS is modified to suppress the unrealistic grid point storms.

(3) The cloud cover calculation is modified which might produce too much low cloud with the modified SC scheme.

(4) The triggering condition is represented by the inclusion of the effects of environmental humidity in the sub-cloud layer and with an upper limit of convective inhibition, which intends to produce more convection in large-scale convergent regions but less convection in large-scale subsidence regions.

The NSAS scheme was implemented to TRAMS by Xu et al. ^[4], and the sensitivities results showed that the NSAS scheme could effectively improve TRAMS performances by using the horizontal exchange of momentum to replace vertical wind shear for parameterization of pressure gradient force. Their results also showed that improvements are made by implementing the exchange coefficient of horizontal momentum that varies with height.

The land-surface model used in the early version of GRAPES included the SLAB scheme (Blackadar ^[21]) and the NOAH scheme (Koren et al. ^[22]). As discussed in Chen et al. ^[1], the SLAB scheme in the early version of GRAPES showed a systematic bias on the simulation of soil flux, which might be caused by the non-inclusion of the changing effects of the soil moisture and diurnal variation of soil temperature in the thin soil layer. For the NOAH scheme, though it included the effects of soil moisture and vegetation canopy, the evaporation was too larger than observation, and the skin sea surface temperature (SST) forecast was not included. To improve the applicability of the land surface model, a simplified model for the land surface scheme was developed with the inclusion of land surface forecasting, SST forecasting, and soil moisture forecasting.

The surface heat balance equation was used to predict soil temperature by taking into account the bottom evolution. The forced-recovery method was used to forecast the soil temperature. The forced-recovery method basically grasps the characteristics of diurnal and annual variations of surface and deep soil temperatures and is still used in many land surface models. A simple "bucket" model is used to predict the surface soil moisture. In addition to the soil temperature and the soil moisture forecasting module, the SMS scheme also included the SST forecasting scheme proposed by Brunke et al. ^[33] According to the evolution of SST, the sensible heat flux and latent heat flux could be calculated accordingly.

Off-line experiments showed that compared with SLAB and NOAH schemes, the sensible heat flux and latent heat flux predicted by SMS schemes were stable and more consistent with observations. Operational experiment results showed that sea and land surface parameterization provided stable and effective predictions of land surface temperature, soil moisture, and SST, which led to an overall improvement of typhoon forecasting.

Although the parameterization of GWD and MBD greatly alleviates the wind bias, the bias still exists in the lower troposphere, in particular over complex terrain and in the PBL. High wind biases are a common phenomenon in most of the models, especially during the simulations over mountains and valleys (Cheng and Steenburgh ^[23]; Skamarock et al. ^[24]; Lorente-Plazas et al. ^[25]). 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 ^[26]). The influences of UTEs include small-scale orography (SSO) effects, which may cause the same order of magnitude as the GWD (Sandu et al. ^[27]).

Zhong et al. ^[6] developed a subgrid orographic parameterization (SOP) scheme which parameterized the UTEs in GRAPES-TMM by adding a sink term in the momentum equations, which was taken as the feedbacks to the momentum tendencies on the first model level in planetary boundary layer parameterization.

The vector of momentum tendency M_V is solved using the matrix A and forcing F_S at the first model level, where A = 1 + σ, and σ represents the sink term. Readers may refer to the introductory part of the SSO parameterization by Zhong et al. ^[6]. Evaluations showed that the surface wind speed bias has been significantly alleviated by adopting the SOP scheme (Fig. 6), in addition to the reduction of the wind bias in the lower troposphere.

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

In the operational evaluation system, the mean track error (units: km) and mean intensity error (hPa) are used in TRAMS, and the traditional evaluation method has been generally used in MARS, including root-mean-square error (RMSE), threat score (TS) for precipitation forecast and bias between forecasts and observations. The mathematical calculation equations are as follows:

where F_C is forecast, O_B is the observation and N is the number of stations in the verification region. OV_bias and UN_bias are the average biases by overestimation and underestimation, respectively.

Shown in Fig. 7 is the annual mean typhoon track error of TRAMS. It can be seen that TRAMS generally presented an annually decreasing trend of typhoon track error. The mean typhoon track error by 24-h, 48-h, and 72-h in 2018 was 77 km, 119 km, and 198 km, respectively. The typhoon intensity error showed a similar decreasing trend as that of the track error. Fig. 8 gives the comparisons of one-year operational evaluations in 2015 by MARS-v1.0 did, and MARS-v2.0. Fig. 8 shows the mean absolute error (AE) of geo-potential height and temperature at a different level in 2015. It can be seen that MARS-v2.0 exhibited much better simulations than v1.0, especially for the geo-potential height under the troposphere (500 hPa).

Figure 7.  Annual mean typhoon track error of TRAMS by 24-h, 48-h and 72-h forecasting.

Figure 8.  Comparisons of mean absolute error (AE) of 48-h simulated geo-potential height and temperature at different levels in 2015 by MARS-v1.0 and MARS-v2.0.

Figure 9 gives the diurnal variations of mean bias of surface temperature and the corresponding number of stations over LG regions. As discussed in Zhong et al. ^[28], MARS showed general underestimation of surface temperature at about 2℃. Both simulations exhibited large sudden changes of the biases during the early morning (0800-1000 LST), which showed a sudden decrease of underestimation of surface temperature by fewer stations and an increase of overestimation of surface temperature by more stations. For the simulations by init-00, the model showed growing number of stations with underestimated surface temperature from morning to late evening, as well as an overwhelmed number of stations with underestimated surface temperature during the nighttime. The underestimation and overestimation reached approximately the same amount in the morning (1000-1200 LST), which might be attributed to the heating effects of short radiation. For the simulations by init-12, the model showed similar characteristics of surface temperature simulation as init-00. However, the underestimated distribution of surface temperature by init-12 was broadly wider than init-00.

Figure 9.  Diurnal variations of mean bias and stations between observation and simulation by underestimation (st_under) and overestimation (st_over) from April to August, 2018: (a) init-00; (b) init-12.

The integrated NWP system at the Guangzhou Regional Meteorological Center plays an important role in the daily weather forecasting over southern China. In particular, the TRAMS model provides a valuable reference for typhoon forecasting. In general, it could predict the typhoon-genesis from 30 hr to 96 hr in advance (Table 1), especially for the typhoon-genesis over South China Sea. Besides, it showed a generally small typhoon-track error and typhoon-intensity error. The 24hr typhoon-track error was between 54 km and 99 km, and the typhoon-intensity error was between 2 hPa and 13 hPa.

Figure 10 demonstrates the capabilities of TRAMS for the typhoon-genesis of two cases over west of the Pacific Ocean. The typhoon "KAJIKI" formed over the South China Sea on 1200 UTC 31 August and the typhoon "LINGLING"formed over the east of the Philippines on 1200 UTC September in 2019. It can be seen that TRAMS could predict the formation of the typhoon "KAJIKI"over the South China Sea at 60 hours in advance and the typhoon "LINGLING"over the east of the Philippines at 96 hours in advance. Both the TRAMS and MARS captured the large-scale distribution and intensity of the outer rainfall over south China by typhoon "KAJIKI"(Fig. 11). However, the model exhibits low predictability of small-scale precipitation, especially that in the warm sector ^[29], which might be attributed to the inaccurate representation of the surface variables (e. g., surface temperature and wind ^{[3, 6, 30]}) and imperfect descriptions of the model physics ^[29], especially for the parameterization effects of the PBL scheme.

Figure 10.  The (a) 60 hr and (b) 96 hr forecasting of sea-level pressure by TRAMS from 1200 UTC 29 August 2019.

Figure 11.  The forecasting of 24 hr accumulative precipitation on 0000 UTC 02 September 2019 by (a) 36 hr forecasting of TRAMS, (b) observation, (c) 36 hr and (d) 48 hr forecasting by MARS.

The present study reviewed the developments of the operational model system based on GRAPES-TMM in GRMC which include the typhoon model of TRAMS, the mesoscale model of MARS, and the rapid update cycling forecasting system. It also reviewed the progress on model dynamics and physics, and evaluated the model with reference to several cases.

In terms of model dynamics, a 3D reference atmosphere scheme has been developed, in which the reference atmosphere is not only a function of height but also a function of longitude and latitude. Moreover, to improve the accuracy of the model, a new scheme of nonlinear terms by fractional steps has been introduced, which first gets the pre-solution of variables before solving the Helmholtz equation, then calculates the nonlinear term and physical feedbacks as an increment of Helmholtz's right-hand term. Besides, based on the quasi-monotone semi-Lagrangian (SL) scheme, the nonlinear constraints used in the linear constraint SL scheme has been developed. Moreover, to explore the effects of the coupling between dynamics and physics on model accuracy, a coupling scheme has been developed which include the feedbacks of temperature and water vapor to atmospheric pressure, which are used to be calculated in the Helmholtz equation of the model as the implicit solution. In general, each improvement of model dynamics contributes to the higher accuracy of the model.

With regard to model physics, parameterization schemes including the subgrid orographic parameterization, a simplified model for the land-surface scheme, a modified cumulus parameterization scheme, and a parameterization scheme of gravity wave drag induced by sub-grid orography have been developed. The operational system has developed other modifications and improvements such as the improvement of the fractional roughness in the PBL scheme, and the adoption of the RRTMG radiation scheme. With the development of observational bases over southern China, the related characteristics and triggering mechanism as well as the development and improvement of the physic package of the model system are further explored.

In terms of model evaluations, both TRAMS and MARS mainly use the traditional evaluation method, including root-mean-square error, threat score for precipitation forecast, and bias between forecasts and observations. Evaluations show that TRAMS presents an annually decreasing trend of typhoon track error, and MARS-v2.0 exhibits much better simulations than v1.0 did. In general, GRAPES-TMM exhibits gradually improved model performances. Note that the traditional evaluation method cannot reflect all the details of the model; more evaluation techniques thus should be included in the future.

It should be noted that the model also exhibits some typical systematic biases such as overestimated wind speed in the troposphere and underestimated nocturnal surface temperature. The simulated time of precipitation onset, peak and subsiding is earlier than the actual time, and the moving speed of convection system is faster than the actual speed. In particular, the model shows some limitations on simulation of the torrential rains in the warm sector over complex terrains. As the higher-resolution upgrading of the NWP modeling system, more model techniques should be further developed. For example, the model dynamics should be calculated with higher accuracy, better vertical stratification and hybrid vertical coordinates. Therefore, the model performance should be further evaluated and improved with more effective descriptions of model dynamics and physic processes.