Climate Modeling with Spectral Elements

1. Climate Modeling with Spectral Elements F. Baer University of MarylandCollege Park, MD 20742 USA Collaborators: H. Wang, J. J. Tribbia and A. Fournier, NCAR M. Taylor, Sandia NOAA NCEP Seminar 24 April 2007

2. Some Current Issues Associated with Climate Modeling • Climate model representation and methodology for computation • Coupling of related models • Parameterization of forces • Use of computing hardware • Speed of computations • Time and space scales required, including regional climate modeling.

3. Primary Goals for this Project and Rationale for using Spectral Elements Improve predictions for weather and climate, global as well as regional, since variable resolution in a model should help in understanding the evolution of the global system. Develop a global atmospheric model methodology which allows for seamless integration of the planetary scales concurrently with regional scales; smaller scale effects often grow on shorter time scales and must be resolved. Optimize calculation time using the features of MPPs to exploit the need for ensemble predictions; Can incorporate the best features of existing state-of-the-art models, in particular the forcing packages.

4. Spectral Elements-Ocean See Haidvogel, et al., Atmosphere-Ocean, 35 (1997) 505-531.

5. The Model for the atmosphere • We have developed a Spectral Element Atmospheric Model (SEAM), a global model with a non-uniform grid offering substantial flexibility and advantages in its features: • It uses geometric properties of finite element methods; • It seamlessly incorporates local mesh refinement with regional detail; • High resolution regions can be applied anywhere on the globe; • It takes advantage of parallel processing; • It maintains the accuracy of spectral models; • It is computationally efficient; • It has no pole problem.

6. Tile a spherical surface with an arbitrary number and size of rectangular elements; Inscribe a polyhedron with rectangular faces inside the sphere, Map the surface of the polyhedron to the surface of the sphere with a gnomonic projection, Use the cube (most elementary polyhedron), Subdivide each of the six cube faces as desired. Use Local Mesh Refinement (LMR) as desired. Define basis functions in each element; Define test functions over the sphere; The Model Domain and Methodology

7. Cube Subdivision 1 Atmosphere Spectral Elements (uniform rectangles)

8. Example of the SEAM Method Shallow water equations Generate an integral form of the equations • Multiply by a global test function ()and integrate over the entire domain (spherical surface). (Note similarity to finite elements.)

9. Spectral element discretization • Decompose sphere into rectangular regions as noted; • Within each rectangle, estimate integral equations by Gauss-Lobatto Quadrature. The Gauss- Lobatto quadrature points on a local rectangle: (8x8 grid is no fluke)

10. Spectral element discretization Within each rectangle, use Legendre cardinal functions for basis functions, , in each direction;

11. Global Test Functions • Simple combinations of Legendre cardinal functions; • One global test function for each grid point. Element boundary points Element interior points • Using the basis and test functions, perform the summation. • Impose continuity of physical functions at element interfaces. • For future tests include flux continuity also.

12. In Summary Tile the sphere with rectangles of arbitrary number and size; Represent the prediction equations in integral form; Use Gauss-Lobatto quadrature for integration; Use Legendre cardinal functions for the basis functions; Use test functions based on the Legendre cardinal functions. These choices result in an extremely simple finite element method with a diagonal mass matrix. The spectral element method appears to be a very efficient and natural way to achieve a high order finite element discretization, and additionally to simplify regional modeling.

13. Local Mesh Refinement (LMR) • Establish a desired uniform global grid; • Choose a regional domain for resolution expansion; • Use a technique called ‘picture framing’ to connect the expanded grid to the courser grid; • We currently use a factor of 3 for this process; • This process may be compounded to get higher expansion regionally. • The process may also be applied to multiple regions.

14. Example of Picture Framing Expanded grid (3x3) (9 elements corner grid (3 elements) Side grid (7 elements)

15. LMR on the globe 3 x LMR over continental USA

16. 8x8 grid per element ~T42 grid LMR over USA (3x3) ~T63 grid

17. Testing the model • Does the method work and can LMR be effectively applied? • Some model test cases with SEAM: • The shallow water equations form: • Use the standard test suite (Williamson et al., 1992); • Use LMR with topography. • The 3-D PE dynamical core (dycore) form; • Use Held-Suarez forcing (simple Newtonian); • Use aqua-planet forcing; • Couple the dycore to full physics and a land surface model from CAM2 (CAM_SEM) and test effects of LMR and an AMIP run; • Compare results against CAM2 runs and data archives where possible.

18. SWE test case 7: Asymmetric polar vortex SE 150x8x8 SE 468x8x8 NCAR T213 SE 1734x8x8 Height field at 5 days

19. Shallow water test case 7: Asymmetric polar vortex • ICs taken from 500 mb observed data; • Compare spectral element model (SE) to NCAR spectral, CSU twisted icosohedral grid (TWIG), and A-L grid point models. • L2 is normalized error between computed and true height field at end of run. SE - spectral element TWIG - twisted icosahedral grid A-L - Arakawa-Lamb

20. Optimal index for Cardinal functions • Curves show increasing number of elements (M) for fixed spectral degree (N); • Increasing N is more costly than increasing M; • N = 8 appears to be the most efficient index. • Select M for desired resolution. Test case 7 Shallow Water Test case 7 M

21. LMR Exp./SWE Test case 5 Grid Global refinement Local refinement 150x8x8 488x8x8 210x8x8 1784x8x8 6936x8x8 310x8x8

22. 488x8x8 LMR Exp. SWE Test case 5 Height field 1784x8x8 Global refinement 6936x8x8 Zoom to Topography 150x8x8 210x8x8 Local refinement 310x8x8

23. Explicit vs Semi-explicit: SWE test case 6 • Errors comparing the explicit and semi-implicit integrations of SEAM. • We also note that the semi-implicit scheme is faster! Explicit Semi-implicit Linf (Taylor) L2 (Thomas and Loft) L1 Days

24. SWE-SEAM / Jupiter A study of decaying turbulence with high resolution using the SWE with SEAM; We use Jupiter dimensions; g = 23 m s-2, radius = 7x104 km One Jupiter day equals 9 Earth hours Equivalent depth = 20000 m Very weak dissipation Dt = 25 Earth seconds T170, T360, T533, T1033 runs on a CRAY T3E with 128 processors T1033 has 60000 elements (8x8) ~ 3000 equatorial pts. Expected result: Rhines scales with 15 jets pole to pole

25. SWE-SEAM / Jupiter • Zonal wind at 700 Jovian days (left column); • Equatorial jet strength in time (right column). • Both demonstrate the effect of resolution. Zonal Wind Equatorial Jet Strength 0 100 T170 200 0 100 T1033 0 200 500 1000 Latitude Jovian days

26. Potential Vorticity on Jupiter T170 T1033 276 Jupiter days

27. Some Demonstrated Features of SEAM Based on the SWE Flexibility: Regional detail and mesh refinement: Advantages on parallel processors: from Jupiter experiment Accuracy: Computational efficiency:

28. The 3-D PE SEAM Dycor Following the SWE tests we created the 3-D PE SEAM dycor (dynamical core); This model was first tested with simple Held-Suarez (H-S) forcing and applied to a polar-jet study; It was next tested to investigate its response to topography, LMR, and various resolutions. These tests indicated that the model was running successfully. The SEAM dycor was then coupled to the NCAR physics package (CAM2) and the land surface model (CLM) and tested by comparing integrations with archived runs made at NCAR with CAM2.

29. Creating CAM_SEM Couple SEAM (dycor) with CAM physics; For initial& boundary data: use a “distance weighted” method to interpolate CAM T42 data to SEAM grids. Couple SEAM with CAM physics (see above) to the CLM; For surface boundary data: use a “distance weighted” method to interpolate the ‘raw’ data to SEAM grids Definitions CAM--Community Atmospheric Model CLM--Community land-surface Model CAM_SEM--Seam dycor with CAM physics and CLM

30. How the couplers work • Tendencies for the dynamical processes are calculated in the dynamical core; i.e., advection of motion variables, etc. • Tendencies for the physics are calculated in the physics package; i.e., radiation heating, convection, etc. • Tendencies for land surface processes are calculated in the CLM model; • When all tendencies are calculated the coupler combines them and extrapolates the variables.

31. Additional Developments Full implementation of - coordinates; Implementation of a semi-Lagrangian scheme for moisture transport; Implementation of an RK4 scheme for time integration with sub-cycling for dynamics: A one-step scheme, consistent with the semi-Lagrange scheme for the spectral element method (Deville et al, 2002). Sub-cycling: use four t/4-steps for dynamics and one t-step for physics.

32. Held-Suarez Experiment Use CAM_SEM dynamical core ( -version) with simple physics and real climatological ICs; Truncate model at T42; Integrate 1200 days with t=450s; Output instantaneous variables every 5 days (interpolate to lat/lon grids); Analyze the output from day 450 to 1200; Compare with output from CAM_EUL run using the same physics and ICs. CAM_EUL--CAM model with Eulerian dycor

33. Held-Suarez Experiment CAM_SEM Pressure levels Height (km) CAM_EUL Zonal Wind (m/s) Temperature

34. Held-Suarez Experiment Height Pressure level CAM_SEM CAM_EUL T- variance (K^2) V-variance ((m/s)^2)

35. Aqua Planet Experiment Use CAM_SEM with CAM physics, but with aqua planet BCs and real climatological ICs; Truncate model at T42; Integrate 1200 days with t=1200s; Output instantaneous variables every 5 days (interpolate to lat/lon grids); Analyze the output from day 450 to 1200; Compare with output from CAM_EUL run using the same physics and ICs. CAM_EUL--CAM model with Eulerian dycor

36. Aqua Planet Experiment Pressure levels Height (km) CAM_SEM CAM_EUL Temperature Zonal Wind (m/s)

37. Aqua Planet Experiment Pressure level Height CAM_SEM CAM_EUL T- variance (K^2) V-variance ((m/s)^2)

38. Full Physics Experiment Use CAM_SEM with CAM physics and CLM; Truncate model at T42; Integrate 1200 days with t=1200s; Output instantaneous variables every 5 days (interpolate to lat/lon grids); Analyze the output from day 450 to 1200; Compare with CAM2.0.1 (NCAR archivedCAM_EULrun with similar conditions), and NCEP reanalysis data. CAM_EUL--CAM model with Eulerian dycor

39. Full Physics ExperimentZonal Wind (m/s) CAM_SEM Pressure levels Height (km) CAM_EUL NCEP REANALYSIS

40. Full Physics ExperimentTemperature (K) CAM_SEM Height (km) Pressure levels CAM_EUL NCEP REANALYSIS

41. LMR Experiment Use CAM_SEM with CAM physics and CLM; Use 3 x3 LMR over the continental USA; ICs and truncation same as Full Physics Exp; Integrate 1200 days with t=600s; Output instantaneous variables every 5 days (interpolated to high-resolution lat/lon grids); Analyze the output from day 600 to 1200; Compare with base (non-LMR) simulation.

42. LMR Experiment Height (km) Pressure levels LMR_USA CAM_SEM Temperature Zonal Wind (m/s)

43. LMR Experiment Height Pressure level LMR_USA CAM_SEM T- variance (K^2) V-variance ((m/s)^2)

44. LMR Experiment LMR_USA Temperature - 200 hPa (K) CAM_SEM

45. LMR Experiment LMR_USA Zonal Wind - 200 hPa (m/s) CAM_SEM

46. LMR Experiment LMR_USA Sea Level Pressure (hPa) CAM_SEM

47. LMR Experiment LMR_USA Total Precipitation Rate (mm/day) CAM_SEM

48. LMR Experiment (ZOOM) LMR_USA Total Precipitation Rate (mm/day) CAM_SEM

49. AMIP Experiment • Use complete CAM_SEM with conditions identical to the CAM2 AMIP archived run. • Integrate from 1979 to 1998; • Set the grid as near T42 as possible; • Compare with output from CAM_EUL runs CAM2, CAM3, and some observations that are available as archives.

50. AMIP Run-Sfc. Pressure