Recent Developments in the NRL Spectral Element Atmospheric Model (NSEAM) * Francis X. Giraldo

1. Recent Developments in the NRL Spectral Element Atmospheric Model (NSEAM)* Francis X. Giraldo www.nrlmry.navy.mil/~giraldo/projects/nseam.html *Funded by the Office of Naval Research

2. Mission of NRL-MRY: Develop NWP models for the Navy Computing Resources: High Performance Computing has shifted towards distributed-memory architectures comprised of tens of thousands of processors (e.g., clustered systems). Current fastest system has 32,000 processors (BlueGene by IBM). Exploiting New Architectures: Need numerical discretization methods that decompose the domain naturally into many smaller subdomains Each subdomain can then be solved independently on-processor Motivation

3. State-of-the-Art Global Forecast Models: Based on the spectral transform method (i.e. spherical harmonics) Spectral transform method uses global approximation functions and Nproc scales like O(T) and costs O(T^3) New Model: NSEAM is based on the spectral element method which combines the high-order accuracy of spectral transforms with domain decomposition property of finite elements Uses local approximation functions and Nproc scales like O(T^2) and costs O(N T^2) Motivation (continued)

4. Domain Decomposition NSEAM: Nproc=O(T^2) ST Models: Nproc = O(T)

5. Spectral Transform Method(Legendre Modal Functions)

6. Spectral Element Method(Legendre Nodal Functions)

7. To construct a model the following things are required: Governing set of Partial Differential Equations along with initial and boundary conditions (PDE) A Coordinate System (CS) A Grid Generator in order to place solution values within the domain (GG) A Spatial Discretization method in order to approximate derivatives (SD) A Time-Integrator in order to advance the solution in time (TI) A Very Large Computer (VLC) What We Need for a New Model

8. Hydrostatic primitive equations (PDE) Equations formulated in Cartesian and sigma coordinates (CS) Spectral Element method in horizontal (SD) Finite Volume method in vertical (SD) 2nd order semi-implicit in time (TI) GMRES with block Jacobi Preconditioning (TI) Message-Passing Interface (VLC) Can use all types of grids including quadrilateral and triangular adaptive grids (GG) Description of NSEAM

9. Hydrostatic Primitive Equations (Mass) (Momentum) (Energy) (State) (Tracers such as Moisture/Salinity)

10. Why Cartesian Coordinates? • Primitive Spherical Form: • Pole problem can be avoided by: • Using Gauss quadrature (as in spherical harmonics) • Creating a hole in the domain (as in gridpoint models) • Using Cartesian coordinates (NSEAM)

11. SE Method • Primitive Equations: • Write Primitive Equations as: • Weak Problem Statement: Find • such that • where

12. SE Method: Basis Functions • Polynomial space for quads and tris • Polynomial Approximation on the Quadrilateral: • With polynomial space • Polynomial Approximation on the Triangle: • With polynomial space

13. Accuracy of SE Method (Shallow Water on Sphere)

14. Grids Icosahedral Lat-Lon Hexahedral Telescoping Thin/Reduced Lat-Lon Adaptive Icosahedral

15. Time-Integrators • For the system: • There are 2 classes of possible time-integrators: • Eulerian - fixed reference frame • Lagrangian - moving reference frame • There are 2 subclasses of these methods: • Explicit – strict time-step restriction • Semi-Implicit – lenient time-step restriction • Fully Implicit – no time-step restriction

16. Explicit Time-Integrator • Start with: • Discretize by 2nd order TI • Everything on the RHS is known. Very fast to solve because only requires the evaluation of a vector of length 3*Np*Nlev (detail of this approach can be found in MWR 2004 paper)

17. Fully-Implicit Time-Integrator • Start with: • Discretize by 2nd order TI • A global sparse NONLINEAR system must be solved. This can be quite expensive! The size of the system is (3*Np*Nlev) x (3*Np*Nlev)

18. Semi-Implicit Time-Integrator • Start with: • Rewrite as: • Discretize by 2nd order TI • A global sparse LINEAR system must be solved. Via a vertical mode decomposition we can do this by solving: Nlev systems of Np x Np sparse linear systems. This is very fast! (details of this approach can be found in QJRMS 2005 paper).

19. Performance for T239 L30 (DT=300 seconds on IBM SP4) (SE Model) (Spectral Model)

20. Performance of NSEAM (216 Processors on IBM SP4) T239 L30 T580 L60 T159 L24

21. Rossby-Haurwitz Wave 4(T=5 days) NOGAPS T80 L24

22. Rossby-Haurwitz Wave 4(T=5 days) NSEAM T80 L24

23. Held-Suarez Test Case(T=1200 days) ECMWF T106 L30

24. Held-Suarez Test Case(T=1200 days) NSEAM T64 L20 Temperature (K) Zonal Wind (m/s)

25. Jablonowski-Williamson Balanced Initial State(Surface Pressure Error for NSEAM T185 L26) Comparison of the explicit (dashed line) and semi-implicit (solid line) NSEAM. The semi-implicit uses a time-step 8 times larger than the explicit. The difference between the two are virtually indistinguishable which is what one would hope to see.

26. Initially-Balanced Still Atmosphere(NSEAM T239 L30 with terrain after 30 days) -90 Meridian 0 Meridian The terrain field for NSEAM T239 from two different views of the globe. In this test an initially-balanced still atmosphere is integrated for 30 days. The atmosphere remained perfectly balanced for the length of the integration which is exactly what should happen.

27. Held-Suarez Test Case(Surface Temp for NSEAM T185 L24 after 30 days) No Terrain Terrain The viewpoint of these plots is from the –120 meridian. The temperature contours look similar for both but the effects of the Rockies, Andes, and Antarctica are clearly visible on the right panel.

28. Jablonowski-WilliamsonBaroclinic Instability(Surface Values at Day 9 for NSEAM T185 L26) Pressure Temperature U Velocity V Velocity

29. Jablonowski-Williamson Baroclinic Instability

30. Jablonowski-Williamson Baroclinic Instability(Surface Values for NSEAM T185 L26 during 0-30 days) Pressure Temperature

31. Future Work • Further testing of NSEAM with terrain • Physical parameterization (sub-grid scale processes) • Non-hydrostatic equations (to better handle dx < 10 km) • Non-hydrostatic Coastal Ocean Model (underway)