a. Weak constraint

Variational assimilation, such as the three-dimensional variational data assimilation (3DVAR) in GSI, attempts to find the solution that best fits observations within a time window to a background (Derber 1989) by minimizing a cost function defined as

View Expanded

wherex is the analysis increment, 𝗕 is the background error covariance,y is the residual observation (difference between the observationy_o and background value 𝗛[x_b]), 𝗥 is the observation error covariance, and 𝗛 represents the observation operators. Note that an additional generic constraint,J_c, has been included. ThisJ_c penalty term can take many forms, but can be designed such that it measures the amount of noise or imbalance generated by the other two terms. Previous studies have incorporated NLNMI operators (Courtier and Talagrand 1990;Thépaut and Courtier 1991) or digital filters (Gauthier and Thépaut 2001;Wee and Kuo 2004) in the design of this term.

A principle problem with various forms of a weak dynamic constraintJ_c is that the minimization problem becomes highly ill-conditioned when the weights forJ_c are increased to a level that produces reasonable balance. Additionally, prescribing the correct weighting of a covariance matrix (analogous to 𝗕 for the background term) for the constraint term can be quite difficult. As a result of these difficulties, the fit of observations to the analysis can become unacceptably degraded compared to the result withoutJ_c.

b. Review of normal-mode initialization

Before describing the idea of the TLNMC as an alternative to a weak-constraintJ_c, a brief review of normal-mode initialization will be presented. The underlying principle of linear and nonlinear normal-mode initialization is that an atmospheric state can be divided into “slow” and “fast” parts, using the normal modes of oscillation of the adiabatic equations of motion linearized about a state of rest. The slow modes correspond closely to traditional geostrophic motion and the fast modes to eastward- and westward-propagating gravity waves. These slow and fast modes were originally used for linear normal-mode initialization (LNMI) byWilliamson and Dickinson (1972) to project the initial atmospheric state onto the gravity modes. These were then subtracted from the original state, supposedly to remove gravity wave oscillations from the forecast. It was later shown byMachenhauer (1977) that LNMI only provided a small reduction in gravity wave noise in the forecast.

Instead of zeroing the gravity mode component of the state,Machenhauer (1977) [and independentlyBaer and Tribbia (1977)] were the first to observe, from empirical and analytic methods, that a much more effective balance could be achieved by projecting the tendencies of the full nonlinear model to gravity modes, and then deriving a correction to the state, which makes the gravity mode tendencies small. In the case of the Machenhauer scheme, a simple Picard iterative technique was used to obtain a solution. This is what makes the procedure nonlinear—thus, NLNMI. This original work was highly successful and resulted in great activity and many variations and applications of NLNMI during the 1980s. But this all came to an abrupt end with the advent of digital filtering (Huang and Lynch 1993), which was much easier to apply with generally more robust results.

Before continuing with a description of the current application of normal-mode initialization methodology, it is natural to ask why use a complex procedure when the much simpler digital filter is just as effective. There are two main reasons: 1) GSI is currently a 3DVAR system and includes a generalized tendency model. The tendency model is used only to estimate an incremental tendency for a single time level and not to do a full time integration, which would be necessary for the application of a digital filter. 2) By projecting directly onto gravity modes, scale selection is possible—an important consideration, since we are primarily interested in removing high-frequency small-scale oscillations, and also some medium fast oscillations with large vertical and horizontal scale. Small-scale, advected Rossby motions can have high local time frequencies, but these are not changed by altering only gravity wave modes in opposition to a single time filtering with a digital filter. So the scale selection is a complex mixture of small and large space scales and medium to high frequency in time. In addition, the digital filter is strictly time selective, and the choice of time periods that are filtered is constrained by the integration period.

a. Description of TLNMC

For those readers interested in details of normal-mode initialization,Daley (1991) is an excellent resource, with extensive references. Here we introduce the following compact notation, which is sufficient to communicate the principle behind TLNMC:

• x_u: unbalanced control analysis increment state vector of lengthn;

• x_c: TLNMC balanced analysis increment state vector of lengthn;

• y: observation innovation vector of lengthp;

• 𝗧:n ×n matrix representing the tangent-linear tendency model mentioned insection 2;

• 𝗙:m ×n matrix, which projects ann-dimensional state vector tom gravity modes (m <n);

• 𝗗:n ×m Machenhauer-style correction matrix used to reduce gravity mode tendencies;

• 𝗖 = 𝗜 − 𝗗𝗙𝗧:n ×n operator equivalent to a single iteration of Machenhauer-type algorithm.

In terms of the above notation, adding TLNMC yields the following modification of Eq.(2.1):

View Expanded

wherex_c = 𝗖x_u is the TLNMC balanced analysis increment andJ_c has been dropped.¹ Note that only the analysis increment is balanced via the TLNMC.

An alternative form for(3.1), expressed entirely inx_c is

View Expanded

where it is assumed that 𝗖 is formally invertible (like 𝗕, the inverse of 𝗖 is never required in the actual minimization algorithm used by GSI). Then the background error covariance for the TLNMC balanced incrementx_c is 𝗖𝗕𝗖^T. Thus, while 𝗕 as currently defined in GSI is fixed and independent of the background state, 𝗖𝗕𝗖^T is a flow-dependent covariance. The addition of a weak-constraintJ_c also modifies the effective background error such that it is also flow dependent. However, the TLNMC operator 𝗖 restricts the unbalanced vectorx_u as a contraction to the observation term. With this contraction, it might be expected that the condition of(3.2) is actually improved over(2.1) with zero weight given to theJ_c term (see theappendix for a more extensive discussion of these points).

It is possible to apply 𝗖 more than once, which would make it truly like the Machenhauer NLNMI, except that the procedure is linear because of the definition of 𝗖. However, this adds significantly to the cost with no apparent benefit over applying 𝗖 just once [see alsoFillion et al. 2007, p. 137, in the context of a limited-area four-dimensional variational data assimilation (4DVAR) analysis system].

b. A practical implementation of TLNMC in GSI

As previously mentioned, to compute the incremental tendencies necessary for the implementation of a NLNMI procedure, a tendency model has been built within the GSI framework. This tendency model is actually a tangent-linear version of a hydrostatic, adiabatic, generalized vertical coordinate, primitive equation model following the formulation inJuang (2005). This model (and the corresponding adjoint) was developed since a simplified, tangent-linear version of the currently operational GFS was not readily available. Additionally, including a generalized tendency model allows for the extension of the TLNMC to non-GFS applications, such as a regional model and domain. The version of the tendency model used for the purposes of a dynamic constraint in GSI is rather simplistic and does not include parameterized physics.

The current version of GSI still uses a library of spherical harmonic transforms to convert the input background state from spectral coefficients to grid variables and the inverse for output of the analysis variables. It was decided that the easiest way to implement TLNMC was to follow precisely the implementation of implicit normal-mode initialization for spectral models (Temperton 1989). This version of NLNMI does not require explicit computation of the horizontal normal modes. Implicit NMI methods, which do not require explicitly defined normal modes, were originally developed for initialization of regional models (Bourke and McGregor 1983;Briere 1982), for which normal modes are hard to define. However, the assumptions made are somewhat more restrictive compared to direct use of normal modes.

The only modification to the spectral scheme presented inTemperton (1989) is to replace the full nonlinear model with a tangent-linear system that is applied to the analysis increment in GSI. It was only necessary to check the algebra in the published paper and then implement it into the GSI code. Thus, the components 𝗙 and 𝗗 of the TLNMC operator 𝗖 were constructed directly from the prescriptions in Temperton.

To render treatment of the 3D system more tractable, the fairly standard treatment of decomposing the solution into vertical modes is employed (Daley 1991;Errico 1991). One problem thus created is that, at each horizontal location, one value of surface pressure andn model levels of temperature are combined inton values of a pseudo geopotential and therefore, when increments of the originaln + 1 values need to be determined from those of the combined field, some closure condition is necessary. This same problem was encountered in the background error definition of the original operational SSI analysis (Parrish and Derber 1992), where the analysis increment of vorticity was statistically coupled to a mass variable, which then needed to be converted to increments of temperature and surface pressure. A constrained variational method was used in the first operational implementation (see appendix A inParrish and Derber 1992), but was later replaced with a direct statistical coupling to temperature and surface pressure because the effective background error variance of temperature was unreasonably large (Derber and Bouttier 1999).

In the current application, a closure scheme that gives a unique, well-behaved solution is used. The method is the same as that described inTemperton and Williamson (1981). This is a dynamically motivated closure consistent with the other assumptions invoked for defining normal modes. An alternative variational approach introduced byDaley (1979) to reduce short-spaced variations (2Δz) in noise in the vertical, which may sometimes occur, was unnecessary in our implementation.

Using vertical modes transforms the 3D problem into a series of 2D shallow-water problems, each with an associated vertical structure function and gravity wave phase speed. The first five vertical modes for a typical case are shown inFig. 1, and the corresponding gravity wave phase speeds are presented inTable 1. These five structures are all approximately barotropic in the troposphere, as a consequence of the high vertical resolution and model top.

In the current formulation of TLNMC in GSI, there are several adjustable parameters. These are 1) NVMODES_KEEP, the number of vertical modes to include, 2) NTLNMC, the number of times to apply the TLNMC operator 𝗖, 3) PERIOD_MAX, the maximum period cutoff in hours for the gravity wave modes to include within each vertical mode system, and 4) PERIOD_WIDTH, also in hours, which specifies the width of a transition zone. This transition zone assures that the horizontal modes are not abruptly cut off, but smoothly vary from used to not used. Note that we do not really have normal modes, just spectral coefficients. A period is assigned to each spectral coefficient based on its two-dimensional wavenumber and vertical mode phase speed. A mask is constructed that is a function of this period and goes smoothly from 1 to 0 as the period increases from PERIOD_MAX − PERIOD_WIDTH to PERIOD_MAX + PERIOD_WITDH. The spectral coefficients of tendencies are first multiplied by this mask, before applying the operator 𝗗. The current settings for operational implementation of TLNMC are NVMODES_KEEP = 8, NTLNMC = 1, PERIOD_MAX = 6 h, and PERIOD_WIDTH = 1.5 h. In theTemperton (1989) formulation, there are three approximations tested, of which system “B” was chosen for the TLNMC application. This scheme was chosen as it is closest to the explicit normal-mode formulation.

a. Single observation tests

The assimilation of a single observation is often used to assess the impact of the assignment of observation error and background error weights, as well as their covariance structures. Because of their simplicity and easy interpretation, such experiments can also be useful for diagnosing potential problems as well as in evaluating analysis changes.Figure 2a shows the resultant 500-hPa temperature analysis increment, for a single 500-hPa temperature observation at 45°N, 90°W with a 1-K innovation and 1-K observation error assimilated with no dynamic constraint. The resultant 500-hPa temperature increment run with the TLNMC turned on (Fig. 2b) is qualitatively and quantitatively similar to the control with no constraint (Fig. 2a). However, the temperature increment is forced to be less circular than the control, increasing the north–south temperature gradient (Fig. 2b). This is also evident when looking at the difference between the TLNMC and control run (Fig. 2c), where a tripole pattern is evident in the difference field. The magnitude of the impact of including the TLNMC appears to be small however, on the order of 5% or less.

For an analysis with no constraint, the assimilation of a single temperature observation results in a wind increment consistent with the definition of the multivariate background error. In the case of the global GSI, the relationship between temperature and wind comes from statistically derived parameters using the so-called NMC method (Wu et al. 2002). A vertical cross section taken at 90°W (Fig. 3a) shows the zonal wind increment determined by the single temperature observation experiment. Through the variable definition and corresponding background error structure functions, a warm temperature observation has resulted in an anticyclonic circulation, with a maximum found near 300 hPa, above the location of the original temperature observation. The same experiment run with the TLNMC turned on increases this anticyclonic circulation substantially (Fig. 3b). The difference between the zonal wind increment for the TLNMC and control runs is as large as 35% at 300 hPa (Fig. 3c).

The fact that the TLNMC increases the induced circulation from a single temperature observation suggests that perhaps the statistically derived background error structure functions are deficient or lacking proper flow dependence. For the midtropospheric midlatitudes, it would be expected that a mass–wind relationship for an analysis increment could be well described as nearly geostrophic, consistent with a thermal wind balance. Though the background error definitions seem to do a reasonable job in capturing a majority of this relationship, the TLNMC appears to clean up the deficiencies in the multivariate background error definition. By first converting a temperature increment to geopotential height increment, one can derive an ageostrophic wind increment by subtracting off the implied geostrophic component from the actual wind increment. For the case without the TLNMC, the single temperature observation resulted in a fairly strong ageostrophic zonal wind increment, particularly above the height of the observation (Fig. 4a). However, the TLNMC seems to be doing what it was designed to do by substantially reducing the unbalanced part of the ageostrophic component of the increment (Figs. 4b,c).

b. Single analysis diagnosis

Further evaluation of the impact of the TLNMC is done by performing a single analysis, utilizing the full suite of observational data, but starting from the same background field. A single analysis was performed for the 1200 UTC GDAS cycle on 9 October 2007, assimilating all operational observations, including conventional and satellite brightness temperature data. The background used for the control and TLNMC assimilation experiments was the operational 6-h forecast from the previous 0600 UTC GDAS cycle. Similar to the single observation experiment, the resultant analysis increments for a control with no constraint and an analysis that utilizes the TLNMC are qualitatively similar (Fig. 5).

The differences that result from using the TLNMC on a single analysis can be quite large however. For the 500-hPa zonal wind increment for this case, the differences are sometimes upward of 50% of the increment (Fig. 6a). Though not shown, this is consistent with other variables and at other levels. Often times, the largest impact of the TLNMC can be found in the extratropics, particularly in regions of strong background flow or regions that contain strong horizontal gradients (Fig. 6). This has important implications, in that the use of such a constraint has in itself the ability to impose some notion of flow dependence. This has also been demonstrated by performing a series of single observation experiments, where the single observation is moved into different parts of the atmosphere where the background has differing characteristics (not shown).

Despite the differences resulting from the TLNMC appearing to be quite large, it is worth noting that the synoptic-scale information content derived from the observations remains intact. This is demonstrated by comparing the midtropospheric meridional wind tendencies in the analyses with and without the TLNMC. The differences in these tendencies are small in magnitude with no discernable systematic pattern evident (Fig. 7). The largest changes to the analysis increment seem to occur in the correction of the barotropic component, which allows for the analysis to fit the observations without removing their synoptic-scale signature, while improving the incremental balance and removing noise.

The TLNMC was designed to remove noise from the analysis while improving the mass–wind relationship where appropriate. One symptom of some of the problems when running the GSI analysis with no constraints is a substantial increase in the amplitude of the zonal mean surface pressure tendency over that of the background (Fig. 8). However, by running the same analysis but with the TLNMC and utilizing the same background as the control assimilation, the problem essentially disappears. The analyzed zonal mean surface pressure tendency from the TLNMC is nearly identical to that derived from the background. This result is consistent with the reduction in the implied (derived) magnitude of the zonal mean RMS of the vertical velocity increment (Fig. 9). This suggests that the analysis from the TLNMC is more in sync with the model dynamics and how the atmosphere behaves, and much less likely to suffer from spinup or spindown [see alsoFillion et al. (2007, their Figs. 6 and 7) for an example of this in the context of a limited-area 4DVAR analysis system where moist physical tangent-linear processes are included].

One of the more exciting elements of the TLNMC is that the improved balance and noise reduction in the analysis increment is achieved without having much of an impact on the minimization itself. The reduction of the norm of the gradient of the cost function for the control analysis without constraint and the TLNMC analysis is nearly identical, with a final reduction of six orders of magnitude [Fig. 10a; see alsoFillion et al. (2007, p. 134), in the context of a limited-area 4DVAR analysis system]. However, this does not all happen without cost, as there does appear to be some impact on the reduction of the cost function (Fig. 10b). The impact on the final reduction of the cost function for this case is approximately 2%. However, this small difference appears to be negligible in the shadows of the positive impact seen from the noise reduction and improved balance and is in fact much smaller than the impact on convergence seen in some weak-constraint experiments.

c. Cycled experiments

Although the single analysis experiments seem to suggest that the quality of analyses is improved through the use of the TLNMC, it is important to evaluate the impact of including such a change on subsequent forecasts. To evaluate the impact of the new constraint on forecast skill, two experiments were run: one without constraint (CTL) and the other with TLNMC turned on (TLNMC). The experiments were designed to mimic NCEP operations, using the operational Global Forecast Model (GFS) and the complete suite of observations as well as the operational data preprocessing codes (as of 1 May 2007). The spectral forecast model was run at T382 truncation with 64 hybrid (sigma pressure) vertical levels, and utilized a digital filter centered about 3 h. The diabatic digital filter was part of all forecast model time integrations, including the shorter 9-h forecasts used for the assimilation as well as the 8-day free forecasts, and was applied in both the CTL and TLNMC experiments. This means that the background for all analyses in both experiments came from model forecasts that had been altered by the digital filter. The GSI analyses were performed on the Gaussian-linear grid that corresponds to the T382 spectral truncation (768 × 384 horizontal grid points) and on the model vertical levels using the operational background error statistics.

Each experiment was fully cycled, creating four analyses and subsequent 9-h forecasts per day respectively. The forecast used as the background for each analysis comes from its own previous cycle. Because the quality control was performed independently for each experiment, it is possible that more data can get into the system as the experiment moves forward in time if the quality of forecasts can somehow be improved. Last, unlike in operations, the full medium-range forecast initialized from the early data cutoff time was only run at 0000 UTC, and not 4 times per simulated day. The experiments were run for 2 months, each starting from the operational GDAS analysis at 0000 UTC 15 November 2006, but the first 2 weeks of results were disregarded for spinup.

It is evident that the use of the TLNMC leads to an improvement of the quality of the short-term forecasts being used as the background for subsequent analysis. The TLNMC experiment shows a consistent reduction in the fit of surface pressure observations to the background, a 6-h forecast, relative to the CTL (Fig. 11). However, the fit of these very same observations to the analysis seems to be slightly less relative to the no-constraint CTL experiment. This slight difference is to be expected as the TLNMC forces the analysis to maintain a sufficiently balanced analysis increment and therefore does not allow the analysis to draw as closely to the individual observations.

The biggest impact on forecast skill, as expected, can be found in the extratropics. In terms of 500-hPa anomaly correlation, the TLNMC is an improvement over the CTL in both hemispheres (Fig. 12). In fact, in the Southern Hemisphere, the use of the TLNMC has increased the forecast skill by as much as one-quarter to one-half day in the medium range. The fact that the gain in the Southern Hemisphere is larger than that seen in the Northern Hemisphere may be seasonal and case dependent. On a forecast by forecast basis, the 5-day 500-hPa anomaly correlation scores from the TLNMC run rarely has dropouts as large as the control, and often significantly outperforms the CTL many days in a row (Fig. 13). In addition to the anomaly correlation, a consistent reduction in the geopotential height root-mean-square error is found (not shown). Results are similar when verification is computed for other levels and variables (not shown).

An improvement in the 3-day precipitation forecasts over North America was also found for the experiment that utilized the TLNMC (Fig. 14). The improvement in skill is small, but consistent, especially for low threshold amounts (i.e., rain versus no rain). The precipitation verification for days 1 and 2 is quite similar, though the differences are even smaller (not shown). However, some degradation in equitable threat score is noted for two of the higher threshold amounts. It is also noteworthy that for day 3, there exists a substantial high bias for the highest threshold amounts (Fig. 14b). However, there are too few cases of the high threshold amounts for this Northern Hemisphere winter period for this to be considered a meaningful result. The TLNMC has been run for other periods (e.g., the Northern Hemisphere summer), and the high bias for heavy precipitation rate threshold amounts does not appear to be an issue.

The analysis differences for the single analysis experiments were quite small in the tropics, when comparing a no-constraint analysis with one that had the TLNMC turned on (Fig. 6b). Although the forecast skill in the tropics appears to be similar between the TLNMC and CTL runs, there does appear to be a small increase in the tropical wind vector RMS error at both 200 (Fig. 15a) and 850 hPa (Fig. 15b). However, this increase is on average only about 1.5% for both pressure levels for a 3-day lead time for this forecast period. Because of the importance of considering diabatic processes in tropical balance, the tangent-linear balance as used here suffers from too weak divergent flow coupled with physical processes. The procedure outlined byFillion et al. (2007) would be a natural extension of the approach taken here, assuming such physical processes could be reliably included into the simplified tendency model.