Search

A study is made of the computational stability of semi-implicit treatments of gravity wave motion suitable for use with two-time-level advection schemes. The analysis is for horizontally uniform reference values of temperature and surface pressure, and for hybrid pressure-based vertical coordinates. Stability requires the use of reference temperatures that are warmer than those that can be used safely with the corresponding three-time-level scheme. The reference surface pressure should also be higher. When stable, the two-time-level scheme is damping, although the largest scales are damped less than by the three-time-level scheme if the latter uses a typical time filtering. The first-order decentered averaging of gravity wave tendencies used in a number of semi-Lagrangian models reduces the need for a relatively warm reference temperature profile but causes a quite substantial damping of otherwise well-represented low-wavenumber modes. The low-wavenumber damping can be avoided by using an alterna...

In this article the implementation of the semi-Lagrangian method in a high-resolution version of the ECMWF forecast model is examined. Novel aspects include the application of the semi-Lagrangian scheme to a global model using the ECMWF hybrid coordinate in the vertical and its use in a baroclinic spectral model in conjunction with a reduced Gaussian grid in the horizontal. The former Eulerian vorticity-divergence formulation is first converted into a momentum-equation formulation that is considerably more economical, thanks in part to the incorporation of Legendre transform efficiencies that were previously demonstrated for the shallow-water equations. The semi-Lagrangian formulation is presented in detail, together with a discussion of computational aspects that are relevant for executing the high-resolution model efficiently on a modestly parallel supercomputer. The impact of formulation changes is assessed via numerical experiments on a set of 12 independent cases. In particular it is shown t...

We compare scalare and vector transform methods for global spectral models of the shallow-water equations. For the scalar transform methods, we demonstrate some economies in the number of Legendre transforms required. It is shown that the vector transform method is algebraically equivalent to the more usual scalar transform methods, and the choice of transform grid is discussed.

In a previous study based on the shallow-water equations, it was shown that nonlinear normal mode initialization (NMI) can be implemented without knowing the normal modes of a model; this implicit form of nonlinear NMI is particularly useful in models for which computing the horizontal normal modes is impracticable. The present paper extends the technique to the multilevel Canadian Operational Finite-Element Regional Model. This paper shows that the method yields well-balanced initial conditions and consistent vertical velocity fields. Forecasts from these initial conditions using a semi-Lagrangian time-integration scheme with relatively large time steps are free from unrealistic high-frequency oscillations.

A three-time-level semi-Lagrangian global spectral model was introduced operationally at the European Centre for Medium-Range Weather Forecasts in 1991. This paper first documents some later refinements to the three-time-level scheme, and then describes its conversion to a two-time-level scheme. Experimental results are presented to show that the two-time-level scheme maintains the accuracy of its three-time-level predecessor, while being considerably more computationally efficient.

The application of variational normal mode initialization to the ECMWF multilevel grid-point model is described. The technique involves minimizing a variational integral of the changes made by the initialization to the analyzed mass and wind fields, suitably weighted with a view to exerting some control over the relative magnitudes of the adjustments to these fields. We show how to construct an appropriate three-dimensional integral, and describe an efficient procedure for solving the minimization problem in the simple case where the weights depend only on latitude. Experimental results are presented to demonstrate that, in comparison with unconstrained initialization, the changes made to the analyzed mass field can be considerably reduced without undue damage to the wind field, and without compromising the benefits of initialization in providing a forecast free of high-frequency oscillations. In the case studied here, the use of constrained normal mode initialization has no significant impact on...

The implementation of fast numerical methods on parallel and vector computers is illustrated by describing the development of fast Fourier transform routines for the vector-processing Cray-1 and Cyber 205 machines. Various vectorization methods are presented for FFT's on the Cray-1. By performing a number of transforms in parallel, “super-vector” performance can be achieved. By modifying the algorithms slightly, multiple transforms can be implemented faster on the Cyber 205 (using 64-bit arithmetic on the 2-pipe model) than on the Cray-1, provided that enough transforms (of order 100) can be performed in parallel.

The semi-implicit semi-Lagrangian integration technique enables numerical weather prediction models to be run with much longer timesteps than permitted by a semi-implicit Eulerian scheme. the choice of timestep can then be made on the basis of accuracy rather than stability requirements. to realize the full potential of the technique, it is important to maintain second-order accuracy in time; this has previously been achieved by applying it in the context of a three-time-level integration scheme. In this paper we present a two-time-level version of the technique which yields the same level of accuracy for half the computational effort. Unlike other efficient two-time-level schemes, ours does not rely on operator splitting. We apply this scheme to a variable-resolution barotropic finite-element regional model with a minimum gridlength of 100 km, using timesteps of up to three hours. the results are verified against a control run with uniformly high resolution, and are shown to be of similar accuracy to those of a semi-implicit Eulerian integration with a timestep of 10 minutes.

A “prime factor” Fast Fourier Transform algorithm is described which is self-sorting and computes the transform in place. With a view to implementation on a Cray-1 or Cyber 205, the form of the algorithm is chosen to minimize the number of additions. With an appropriate choice of index mapping in the derivation, we obtain the unexpected result that the required indexing is actually simpler than that for a conventional FFT. The construction of the necessary “rotated” DFT modules is described, and comparisons are presented between the new algorithm and the conventional FFT in terms of operation counts and timings on an IBM 3081; on this machine, the new transform algorithm takes about 60% of the time for the conventional FFT. A Fortran routine for the new algorithm is outlined.

Orthogonal vertical normal modes are found for a multilevel sigma-coordinate model, linearized about a state of rest with a nonisothermal mean temperature profile. Orthogonality permits the partitioning of energy into the vertical modes, and simplifies the application of variational techniques to normal mode initialization in the multilevel case. Improved procedures are described for the vertical transforms required during normal mode initialization. We derive relationships between the time derivatives of the energy computed in vertical normal mode space and in gold point space, and analogous relationships between the changes made by initialization to the vertical normal mode coefficients and the corresponding changes made in grid point space.

It is shown that the self-sorting variants of the mixed-radix FFT algorithm may be specialized to the case of real or conjugate-symmetric input data. In comparison with conventional procedures, savings of around 20% are achieved in terms of operation counts. A multiple real/half-complex transform package on the Cray-1, based on the algorithms described here, achieves a 30% saving in CPU time compared with a package using conventional algorithms. A similar package has also been implemented on the Cyber 205.

We present the firm application of the semi-implicit semi-Lagrangian integration technique to a finite-element barotropic model using the shallow-water equations on a variable-resolution. Two schemes based on this approach, differing only in their treatment of the rotational part of the wind field, are formulated and analyzed. A set of comparative experiments was performed using carefully balanced initial conditions (to eliminate spurious high-frequency oscillations); an Eulerian control integration was run at high resolution on a uniform grid. Both schemes are stable with timesteps at least six times longer than the limiting timestep of the corresponding Eulerian scheme using the same variable-resolution mesh. However, one scheme is consistently more accurate than the other. These results were explained by a theoretical analysis of the stability and accuracy of the schemes. We conclude that the semi-implicit semi-Lagrangian scheme is a promising technique for finite-element models as well as for...

One of the simplest fast direct methods for solving the discrete Poisson equation is adapted for use in very large problems where the right-hand side and solution fields must be stored on disc. The performance of the algorithm described here is compared with that of a technique due to Schumann for a three-dimensional Poisson problem on a (128)3 grid.

It is shown that nonlinear normal mode initialization (NMI) can be implemented without knowing the normal modes of a model. The implicit form of nonlinear NMI is particularly useful for models whose normal modes cannot readily be found; for example, if the underlying linear equations are nonseparable. An implicit nonlinear NMI scheme is formulated for the shallow-water equations on a polar stereographic projection. The linear equations which define the implicit normal modes include most of the beta terms as well as variable Coriolis parameter and map scale factor. Even in this nonseparable case, the equivalence between implicit and conventional nonlinear NMI is shown to be exact. The scheme is implemented in a regional model on a quasi-hemispheric domain, which uses a finite-element discretization on a nonuniform grid. The well-posed lateral boundary conditions of this model lead to consistent boundary conditions for the initialization. Results are presented not only for the implicit form of Mach...

Direct methods for the solution of the discrete Poisson equation over a rectangle are commonly based either on Fourier transforms or on block-cyclic reduction. The relationship between these two approaches is demonstrated explicitly, and used to derive the FACR( l ) algorithm in which the Fourier transform approach is combined with l preliminary steps of cyclic reduction. It is shown that the optimum choice of l leads to an algorithm for which the operation count per mesh point is almost independent of the mesh size. Numerical results concerning timing and round-off error are presented for the N × N Dirichlet problem for various values of N and l . Extensions to more general problems, and to implementation on parallel or vector computers are briefly discussed.

The author shows that prime factor FFT algorithms offer little improvement over conventional FFT algorithms on computers such as the Cray-1 and Cyber 205 where the multiplications can be performed in parallel with the additions. A very modest gain may be obtained by using Good's algorithm (1958) with conventional small-n transforms. Winograd's technique (1978), despite its impressive reduction in the number of multiplications, is likely to be slower than the conventional algorithm, particularly on the Cray-1, where memory transfers will dominate the computation. 13 references.

Some comparisons are presented between various algorithms for solving the discrete Poisson equation over a rectangle with Dirichlet boundary conditions on two opposite sides and either periodic or Dirichlet boundary conditions on the other two sides. The methods considered include those based on FFT in one-dimension, block-cyclic reduction, and variants of the FACR(1) algorithm in which FFT is combined with a preliminary cyclic reduction step. Comparisons are made in terms of detailed operation counts, execution time, and accuracy with respect to round-off error, for a wide range of grid sizes. The importance of careful programming is stressed, and some suggestions are offered for efficient implementation of the algorithms.

A new fast Fourier transform algorithm for real or half-complex (conjugate-symmetric) input data is described. Based on the decomposition of N (the length of the transform) into mutually prime factors, the algorithm performs transforms in-place and without pre- or post-reordering of the data. With large-scale scientific computing in mind, the emphasis is on reducing the number of additions required. Compared with the best available algorithm based on specializing the conventional FFT to the real/half-complex case, the number of multiplications is also reduced by about 50%. On the Cray X-MP, a transform package based on the new algorithm runs up to 20% faster than the previous fastest available routines, besides halving the storage requirements.

Implicit nonlinear normal mode initialization schemes enable nonlinear NMI to be performed in models whose normal modes cannot readily be computed. Such schemes are algebraically equivalent to conventional (“explicit”) NMI based on the normal modes of a set of linearized equations which is slightly different from the usual choice. In this paper we apply implicit NMI to a barotropic spectral model whose normal modes are easily found, thus permitting a direct comparison to be made between conventional and implicit NMI schemes. Both first-order (Machenhauer) and second-order (Tribbia) variants of implicit nonlinear NMI are formulated for a spectral model and compared with their explicit counterparts. Experimental results show that the differences between explicit and implicit nonlinear NMI am insignificant except at the very largest horizontal scales. Besides validating the concept of implicit nonlinear NMI, this study suggests a practical approach to initializing very high resolution spectral model...

Some examples are given of the uses of the FFT algorithm in numerical weather prediction. The algorithm is derived by means of matrix factorization, in such a way that the many possible variants can be obtained in a unified manner. Particular attention is paid to self-sorting variants of the algorithm, but skeleton Fortran routines are included for all the “canonical” forms. The FFT works by decomposing the discrete Fourier transform into a sequence of short transforms. Algorithms for these short transforms are presented and used to derive operation counts for the whole algorithm. The implementation of the FFT on vector computers is described, and in the final section it is demonstrated how savings can be achieved in the case of two-dimensional transforms.

It is shown that implicit normal mode initialization can be combined with a variational technique, in order to control the relative magnitudes of the changes to the analyzed mass and wind fields. Since the initialization procedure is expressed entirely in physical space, the use of locally varying weights in the variational integral becomes more straightforward than in previous efforts to combine variational methods with normal mode initialization. We present details of the application to a finite-element model of the shallow water equations on a stereographic projection. It is demonstrated that the use of variational initialization can change the slowly evolving component of the subsequent forecast, as well as eliminate the unrealistic fast component.

In Part I of this paper we review initialization methods for numerical weather prediction models, leading up to the development of schemes based on the normal modes of the forecast model. We present the derivation of the normal modes of ECMWF's multilevel global grid-point model, and compare the horizontal normal modes with those obtained using alternative finite-difference schemes. The impact of stability-enhancing Fourier filtering procedures on the normal modes is also discussed. Finally in Part I we apply linear normal mode initialization to a nine-level version of the model with 3.75° horizontal resolution. The application of nonlinear normal mode initialization to this model is presented in Part II.

An implicit version of nonlinear normal mode initialization (NMI) is briefly described; this procedure enables nonlinear NMI to be performed even in numerical weather prediction models for which the normal modes cannot readily be found. The development of this implicit scheme also considerably simplifies the application of variational constraints to the initialization procedure.

During the past five years I have had the good fortune to be able to work on two of the world's fastest computers, the Cray-1 (at the European Centre for Medium Range Weather Forecasts) and the Cyber 205 (at the UK Meteorological Office). A comparison between the two machines is clearly of interest not only to potential users but also to a more general audience. In this article I will point out some of the differences between the two machines, and the impact these differences have on the way in which certain problems should be formulated and programmed.