Adiabatic formulation of the ECMWF model

1. Adiabatic formulation of the ECMWF model Agathe Untch e-mail: Agathe.Untch@ecmwf.int (office 11) Adiabatic formulation of the ECMWF model

2. Introduction • Step by step guide through the decisions to be taken / choices to be made when designing the adiabatic formulation of a global Numerical Weather Prediction (NWP) model. • In the process we are “constructing” the dynamical core of the ECMWF operational NWP Model. Adiabatic formulation of the ECMWF model

3. Introduction (cont.) • A numerical model has to be: • stable • accurate • efficient • No compromise possible on stability! • The relative importance given to accuracy versus efficiency depends on what the model is intended for. • For example: • an operational NWP model has to be very efficient to allow the running of all applications (data-assimilation, forecasts, ensemble prediction system) in a tight daily schedule. • a research model might not have to be so efficient but can’t compromise on accuracy. Adiabatic formulation of the ECMWF model

4. Introduction (cont.) • Essential to the performance of any NWP or climate-prediction model are a.) the form of the continuous governing equations (approximated or full Euler equations?) b.) boundary conditions imposed (conservation properties depend on these). c.) the numerical schemes chosen to discretize and integrate the governing equations. Adiabatic formulation of the ECMWF model

5. Euler Equations for a moist atmosphere on a rotating sphere 3D momentum equation Continuity equation Thermodynamic equation Humidity equation Transport equations of various physical/chemical species Equation of state Adiabatic formulation of the ECMWF model

6. specific volume total time derivative Notations: q specific humidity L latent heat Xi mass mixing ratios of physical or chemical species (e.g. aerosols, ozone) g gravity = gravitation g* + centrifugal force Spherical geopotential approximation is made: neglect Earth’s oblateness (~0.3%). => spherical geometry assumed! Adiabatic formulation of the ECMWF model

7. Momentum equations in spherical coordinates : With Euler Equations in spherical coordinates x Adiabatic formulation of the ECMWF model

8. Continuity equation in spherical coordinates With With Thermodynamic equation in spherical coordinates Adiabatic formulation of the ECMWF model

9. Shallow atmosphere approximation 1. Replace r by the mean radius of the earth a and by , where z is height above mean sea level. 3. Neglect all metric terms not involving 4. Neglect the Coriolis terms containing (resulting from ). the horizontal component of a a a z 2. Neglect vertical and horizontal variations in g. a Adiabatic formulation of the ECMWF model

10. Euler Equations in shallow atmosphere approximation is horizontal velocity n is the tracer for the hydrostatic approximation: n=0 => vertical momentum equation = hydrostatic eq. For n=1 we refer to these equations as Non-hydrostatic equations in Shallow atmosphere Approximation (NH-SA) Adiabatic formulation of the ECMWF model

11. Choice of predicted variables Predicted variables: ? or => (Allows enforcement of mass conservation) (Thermodyn. Computations are simpler) Combine continuity, thermodynamic & gas equation to obtain a prognostic equation for p. Form of NH-SA equations more commonly used in meteorology: Adiabatic formulation of the ECMWF model

12. Hydrostatic Approximation (n=0) • Benefits from hydrostatic approximation • Vertical momentum equation becomes a diagnostic relation (=> one prognostic variable (w) less!) • Vertically propagating acoustic waves are eliminated (these are the fastest waves in the atmosphere, causing the biggest stability problems in numerical integrations!) • Drawbacks of hydrostatic approximation • Not valid for short horizontal scales (for mesoscale phenomena) • Short gravity waves are distorted in the hydrostatic pressure field. • Operational version of the ECMWF model is a hydrostatic model. • operational horizontal resolution ~25km (T799), so hydrostatic approximation is (still) OK. Adiabatic formulation of the ECMWF model

13. Hydrostatic shallow atmosphere equations(Hydrostatic Primitive Equations (HPE)) • pis monotonic function ofz and can be used as vertical coordinate. (Eliassen (1949)) Adiabatic formulation of the ECMWF model

14. Potential temperature (isentropic coordinate): • good coordinate where atmosphere is stably stratified • (potential temperature increases monotonic with z). • adiabatic flow stays on isentropic surfaces (2D flow) • good coord. in stratosphere, not very good in troposphere Choice of vertical coordinate Most commonly used vertical coordinates: Height above mean sea level z: - most natural vertical coordinate Pressure p (isobaric coordinate): - has advantages for thermodynamic calculations - makes the continuity equation a diagnostic relation in the hydrostatic system - can be extended for use in non-hydrostatic models Adiabatic formulation of the ECMWF model

15. Choice of vertical coordinate (cont.) Generalized vertical coordinate s: (Kasahara (1974), Staniforth & Wood (2003)) Any variable s which is a monotonic single-valued function of height zcan be used as a vertical coordinate. Coordinate transformation rules (from z to any vertical coordinate s): Adiabatic formulation of the ECMWF model

16. Hydrostatic relation between p and z: => with geopotential Pressure p as vertical coordinate in the hydrostatic system Coordinate transformation rules (from z to p): Adiabatic formulation of the ECMWF model

17. Hydrostatic Primitive Equations with pressure as vertical coordinate Pressure gradient replaced by geopotential gradient (at constant pressure). pressure vertical velocity Continuity eq. is a diagnostic eq. in p-coordinates. • Number of prognostic variables reduced to 3 (horizontal winds & T)! Geopotential computed from hydrostatic equation. Adiabatic formulation of the ECMWF model

18. For every s for which => continuity is diagnostic eq.! Hydrostatic pressure as vertical coordinate for a non-hydrostatic shallow atmosphere model Introduced by Laprise (1991) In the hydrostatic system with pressure as vertical coordinate the continuity equation is a diagnostic equation. The idea is to find a vertical coordinate for the NH system which makes the continuity equation a diagnostic equation. Continuity equation in generalized vertical coordinate s: (Kasahara(1974)) Adiabatic formulation of the ECMWF model

19. Choose and denote with the coordinate s for which i.e. For is the weight of a column of air (of unit area) above a point at height z, i.e. hydrostatic pressure. Hydrostatic pressure as vertical coordinate for a non-hydrostatic shallow atmosphere model (cont.) Adiabatic formulation of the ECMWF model

20. Hydrostatic pressure as vertical coordinate for a non-hydrostatic shallow atmosphere model (cont.) => D3 inZ Adiabatic formulation of the ECMWF model

21. Boundary conditions Governing equations have to be solved subject to boundary conditions. • The lower boundary of the atmosphere (surface of the earth) is • a material boundary (air parcel cannot cross it!) • velocity component perpendicular to surface has to vanish (e.g. at a flat and rigid surface vertical velocity w = 0) Unfortunately, the topography of the earth is far from flat, making it quite tricky to apply the lower boundary condition. Solution: Use a terrain-following vertical coordinate. For example: traditional sigma-coordinate (Phillips, 1957) Adiabatic formulation of the ECMWF model

22. Terrain-following vertical coordinate The easiest case is a flat and rigid boundary where the boundary condition simply is: e.g. at the bottom & at the top We are looking for a vertical coordinate “s” which makes it easy to apply the condition of zero velocity normal to the boundary, even for very complex boundaries like the earth’s topography. at the boundary. Therefore, we look to create a vertical coordinate which makes the upper and lower boundaries “flat”. That is, s is constant following the shape of the boundary (i.e. the boundary is a coordinate surface). Adiabatic formulation of the ECMWF model

23. Terrain-following vertical coordinate(cont.) is a monotonic single-valued function of hydrostatic pressure and also depends on surface pressure in such a way that (Where is the pressure at the top boundary.) this is the traditional sigma-coordinate of Phillips (1957) For The ECMWF model uses a terrain-following vertical coordinate based on hydrostatic pressure. The principle will be explained based on hydrostatic pressure : A simple function that fulfils these conditions is Adiabatic formulation of the ECMWF model

24. Sigma-coordinate (First introduced by Phillips (1957)) Drawback: Influence of topography is felt even in the upper levels far away from the surface. Remedy: Use hybrid sigma-pressure coordinates x Adiabatic formulation of the ECMWF model

25. Hybrid vertical coordinate First introduced by Simmons and Burridge (1981). The difference to the sigma coordinate is in the way the monotonic relation between the new coordinate and hydrostatic pressure is defined: The functions A and B can be quite general and allow to design a hybrid sigma-pressure coordinate where the coordinate surfaces are sigma surfaces near the ground, gradually become more horizontal with increasing distance from the surface and turn into pure pressure surfaces in the stratosphere(B=0). In order that the top and bottom boundaries are coordinate surfaces (=> easy application of boundary condition), A and B have to fulfil: at the surface at the top Adiabatic formulation of the ECMWF model

26. Comparison of sigma-coordinates & hybrid η-coordinates Coordinate surfaces over a hill for sigma-coordinate η-coordinate Adiabatic formulation of the ECMWF model

27. Non-hydrostatic equations in hybrid vertical coordinate prognostic continuity eq. Adiabatic formulation of the ECMWF model

28. Hydrostatic Primitive Equations in hybrid ηvertical coordinate In addition to the geopotential gradient term, a pressure gradient term again! Continuity equation is prognostic again because the (hydrostatic) pressure is not the vertical coordinate anymore. Adiabatic formulation of the ECMWF model

29. Hydrostatic equations of the ECMWF operational model(incorporating moisture) Adiabatic formulation of the ECMWF model

30. virtual temperature gas constant of dry air, gas constant of water vapour specific heat of dry air at constant pressure specific heat of water vapour at constant pressure p-coordinate vertical velocity Notations: q specific humidity Xmass mixing ratio of physical or chemical species (e.g. aerosols, ozone) contributions from physical parametrizations horizontal diffusion terms Adiabatic formulation of the ECMWF model

31. with boundary conditions: B from def. of vert. coord. From the continuity equation we can derive (by vertical integration) the following equations: Needed for the energy-conversion term in the thermodynamic equation Needed for the semi-Lagrangian advection Prognostic equations for surface pressure Adiabatic formulation of the ECMWF model

32. Prognostic equations of the ECMWF hydrostatic model These equations are discretized and integrated in the ECMWF model. Adiabatic formulation of the ECMWF model

33. Discretisation in the ECMWF Model • Space discretisation • In the horizontal: spectral transform method • In the vertical: cubic finite-elements • Time discretisation • Semi-implicit semi-Lagrangian two-time-level scheme Adiabatic formulation of the ECMWF model

34. Horizontal discretisation Options for discretisation are: • in grid-point space only (grid-point model) • finite-difference, finite volume methods • (in spectral space only) • in bothgrid-point and spectral space and transform back and forth between the two spaces (spectral transform method, spectral model) • Gives the best of both worlds: • Non-local operations (e.g. derivatives) are computed in spectral space (analytically) • Local operations (e.g. products of terms) are computed in grid-point space • The price to pay is in the cost of the transformations between the two spaces • in finite-element space(basis functions with finite support) Adiabatic formulation of the ECMWF model

35. Horizontal discretisation (cont.) ECMWF model uses the spectral transform method Representation in spectral space in terms of spherical harmonics: Ideally suited set of basis functions for spherical geometry (eigenfunctions of the Laplace operator). m: zonal wavenumber n: total wavenumber λ= longitude μ= sin(θ) θ: latitude Pnm: Associated Legendre functions of the first kind Adiabatic formulation of the ECMWF model

36. The horizontal spectral representation Triangular truncation (isotropic) Spherical harmonics Fourier functions associated Legendre polynomials FFT (fast Fourier transform) using NF 2N+1 points (linear grid) (3N+1 if quadratic grid) Legendre transform by Gaussian quadrature using NL (2N+1)/2 “Gaussian” latitudes (linear grid) ((3N+1)/2 if quadratic grid) No “fast” algorithm available Triangular truncation: n N m m = -N m = N Adiabatic formulation of the ECMWF model

37. Horizontal discretisation (cont.) Representation in grid-point space is on the reduced Gaussian grid: Gaussian grid: grid of Guassian quadrature points (to facilitate accurate numerical computation of the integrals involved in the Fourier and Legendre transforms) - Gauss-Legendre quadrature in latitude: Grid-points in latitude are the zeros of the Legendre polynomial of order NG Gaussian latitudes NG (2N+1)/2 for the linear grid. NG (3N+1)/2 for the quadratic grid. - Gauss-Fourier quadrature in longitude: Grid-points in longitude are equidistantly spaced (Fourier) points 2N+1 for linear grid 3N+1 for quadratic grid Adiabatic formulation of the ECMWF model

38. The Gaussian grid Reduced grid Full grid About 30% reduction in number of points • Associated Legendre functions are very small near the poles for large m Adiabatic formulation of the ECMWF model

39. T799 T1279 25 km grid-spacing ( 843,490 grid-points) Current operational resolution 16 km grid-spacing (2,140,704 grid-points) Future operational resolution (from end 2009) Adiabatic formulation of the ECMWF model

40. Spectral transform method Grid-point space -semi-Lagrangian advection -physical parametrizations Inverse FFT FFT Fourier Space Fourier Space Spectral space -horizontal gradients -semi-implicit calculations -horizontal diffusion Inverse LT LT FFT: Fast Fourier Transform, LT: Legendre Transform Adiabatic formulation of the ECMWF model

41. Horizontal discretisation (cont.) Advantages of the spectral representation: a.) Horizontal derivatives are computed analytically => pressure-gradient terms are very accurate => no need to stagger variables on the grid b.) Spherical harmonics are eigenfunctions of the the Laplace operator => Solving the Helmholtz equation (arising from the semi-implicit method) is straightforward. => Applying high-order diffusion is very easy. Disadvantage: Computational cost of the Legendre transforms is high and grows faster with increasing horizontal resolution than the cost of the rest of the model. Adiabatic formulation of the ECMWF model

42. Comparison of cost profilesat different horizontal resolutions Adiabatic formulation of the ECMWF model

43. Cost of Legendre transforms T2047 T1279 T799 T511 Adiabatic formulation of the ECMWF model

44. Profile for T2047on IBM p690+ (768 CPUs) Legendre Transforms ~17% of total cost of model Physics ~36% of total cost Adiabatic formulation of the ECMWF model

45. L91 L60 Vertical discretisation Variables are discretized on terrain-following pressure based hybrid η-levels. Vertical resolution of the operational ECMWF model: 91 hybrid η-levels resolving the atmosphere up to 0.01hPa (~80km) (upper mesosphere) Adiabatic formulation of the ECMWF model

46. Vertical discretisation (cont.) Choices: - finite difference methods - finite element methods Operational version of the ECMWF model uses a cubic finite-element (FE) scheme based on cubic B-splines. No staggering of variables, i.e. all variables are held on the same vertical levels. (Good for semi-Lagrangian advection scheme.) Inspection of the governing equations shows that there are only vertical integrals (no derivatives) to be computed (if advection is done with semi-Lagrangian scheme). Adiabatic formulation of the ECMWF model

47. Prognostic equations of the ECMWF hydrostatic model Reminder: slide 32 These equations are discretized and integrated in the ECMWF model. Adiabatic formulation of the ECMWF model

48. with boundary conditions: B from def. of vert. coord. Reminder: slide 31 From the continuity equation we can derive (by vertical integration) the following equations: Needed for the energy-conversion term in the thermodynamic equation Needed for the semi-Lagrangian advection Prognostic equations for surface pressure Adiabatic formulation of the ECMWF model

49. Aji Bji Vertical integration in finite elements can be approximated as Basis sets Applying the Galerkin method with test functions tj => Adiabatic formulation of the ECMWF model

50. Vertical integration in finite elements Including the transformation from grid-point (GP) representation to finite-element representation (FE) and the projection of the result from FE to GP representation one obtains Matrix J depends only on the choice of the basis functions and the level spacing. It does not change during the integration of the model, so it needs to be computed only once during the initialisation phase of the model and stored. Adiabatic formulation of the ECMWF model