340 likes | 486 Vues
Recent developments in Numerical Methods - with a report from the ECMWF annual seminar on ‚ Recent developments in numerical methods for atmosphere and ocean modelling 2013‘. 29 th WGNE meeting 10-14 March 2014, Melbourne M. Baldauf (DWD).
E N D
Recentdevelopments in NumericalMethods -with a reportfromtheECMWFannualseminar on ‚Recentdevelopments in numericalmethodsforatmosphereandoceanmodelling 2013‘ 29th WGNE meeting 10-14 March 2014, Melbourne M. Baldauf (DWD) M. Baldauf (DWD)
ECMWF annualseminar 2013 on ‚Recentdevelopments in numericalmethodsforatmosphereandoceanmodelling‘ 2-5 Sept. 2013, ECMWF, Reading, UK organizer: Nils Wedi (ECMWF) invitedtalks also frommanygloballymodellingcenters http://www.ecmwf.int/newsevents/meetings/annual_seminar/2013/index.html (last annualseminars on numericalmethods: 2004, 1998, 1991, …) horizontal grid: N. Wood, G. Zängl, B. Skamarock, H. Tomita, J. Thuburn vertical discret.: G. Zängl, M. Hortal, M. Baldauf Finite element based: C. Cotter, F. Giraldo, M. Baldauf time integration: P. Benard, R. Klein, S.-J. Lock unstructured meshes: J. Szmelter reduced equation systems: P. Smolarkiewicz, R. Klein Physics-dynamics-coupling: S. Malardel M. Baldauf (DWD)
Classicalrequirementsfor a dynamicalcore: • Stability • Accuracy (upto a certain order; at least 2) • Efficiency • Today, additional requirements: • conservationofcertain variables • „mimetic“ properties (discretizationreproducessomeexactanalyticalproperties) • well-balancing • scalability • efficient on (quite) different computerarchitectures (CPU, GPU,…) M. Baldauf (DWD)
Development and projection of HPC over the years (source: www.top500.org) 1 EFlops 3*106 cores 1 PFlops sum over all 500 estimated ~108 cores ~104 cores #1 1 TFlops factor 1.8 per year = factor 2 per 14 months #500 1 GFlops code scalability is crucial 1993 2003 2013 M. Baldauf (DWD)
Since the clock rate of supercomputers does not increase increase of computer power only by parallelisation scalability of the code is very important Many current global models (still) use a lat-lon grid clustering of grid points around the pole get rid of the strong time step restriction by combination of Semi-Lagrangian-advection and semi-implicit time integration SL for high Courant numbers requires intensive communication detrimental for scalability look for horizontally more uniform grids (lat-lon is not a problem for limited area models) M. Baldauf (DWD)
Horizontal grids review of horizontal grids on the sphere: Staniforth, Thuburn (2012) QJRMS • Whybotheringwith all this horizontal gridstuff? • more uniform grids (seeabove) • gridstaggering • degreesoffreedom (DoFs) per celland/or variable • see also a generaloverviewabout ‚computationalmodes‘ by J. Thuburn M. Baldauf (DWD)
Why grid staggering? Dispersion relation of the linear 1D wave equation Frequency Phase-, group-velocity =0 for 2x waves x/c vph unstaggered negative group velocity vgr kx kx grid tstagg = ½ tunstagg x/c vph staggered vgr kx kx M. Baldauf (DWD)
N. Wood (UK MetOffice) M. Baldauf (DWD)
N. Wood (UK MetOffice) M. Baldauf (DWD)
unfortunately low accuracy of discretizations N. Wood (UK MetOffice) M. Baldauf (DWD)
B. Skamarock (NCAR) M. Baldauf (DWD)
Currentlyusedgridsforoperational models • UK MetOffice: currentlylat-longrid (‚New Dynamics‘, ‚ENDGame‘)GungHo: more uniform grids (exampleseeabove) • NCAR: MPASusesicosahedron dual grid hexagonal C-grid • DWD: currently (GME): icosahedron (however, regulargrid on 10 diamonds)new: ICON usesicosahedronwithtrianglerefinement, triangular C-grid • Canada: currentlylat-lonplans: Yin-Yang grid in 2015 • China (CMA): FD on regularlat-longrid SISLschemenew: (multi-moment) FD; Yin-Yang gridorcubedsphere. • … M. Baldauf (DWD)
B. Skamarock (NCAR) this is still a conformal grid (i.e. no hanging nodes) M. Baldauf (DWD)
alternatively: hide most of the problems with the horizontal grid behind a spectral representation • ECMWF: horizontallyspectral (verticallyFE) • MeteoFrance: currentlyspectralmodels (sharedyn. corewithECMWF)however, visionforFEdiscretizationhorizontally • JMA: GSM usesspectralapproach • Brazil: spectral / lat-lon SI integration • … M. Baldauf (DWD)
The Integrated Forecasting System (IFS) technology applied at ECMWF for the last 30 years … A spectral transform, semi-Lagrangian, semi-implicit (compressible) (non-)hydrostatic model Big obstacle could be removed: Legendre transform is O(N^3), fast Legendre Transform O(N^2 log N) (Wedi et al. (2013), Tygert (2008, 2010) ) Similar statements hold for Aladin (biperiodic Fourier series) N. Wedi (ECMWF)
N. Wedi (ECMWF) M. Baldauf (DWD)
by G. Mozdzynski M. Baldauf (DWD)
Vertical grids • issues: • vertical finite elementrepresentation • proper verticalaverages in verticallystaggeredgridforFDarenecessary M. Baldauf (DWD)
Improvement of the vertical discretization Vertical (Lorenz)-staggering in COSMO (and a few other models): Half levels (w-positions)are defined by a stretching function zk = f(k); Main levels (p‘, T‘-pos.)lie in the middle of two half levels Arithmetic average from half levels to main level: Weighted average from main levels to half level Derivatives always by centered differences (appropriate average used before) The above mentioned staggering requires proper averaging; example … M. Baldauf (DWD)
Truncation error in stretched grids Divergence – grid stretching variant B Divergence with weighted average of u (andv) to the half levels: s·dz dz 1/s ·dz Divergence with only arithmetic average of u (and v) to the half levels: not a consistent discretization if s1 ! M. Baldauf (DWD)
Time integration • Meteo France (P. Benard): 3-time level Semi-implicitLeapfrogbasedintegration+ iterationstep (ICI-scheme) (on spectralmodel) • DWD (ICON): 2-time levelPredictor-Corrector –scheme (HEVI in eachstep)COSMO: 2-timelevelRunge-Kutta (split explicit) • NCAR: 2-timelevel Runge-Kutta (split explicit) forWRFandMPAS • UK MetOffice: currently SI-SL schemeoutlook: IMEX Runge-Kutta underconsideration • JMA: GSM uses SI-SL scheme (forthespectral model) • … M. Baldauf (DWD)
Advectionschemes • Fortraceradvectionbeyondsufficientaccuracy, twoingredientsareseenasimportant: • conservation finite volumediscretization • consistencywithcontinuityequation(=theadvectionequationforq shouldreproducethoseof if q=1)easy, ifthe same advectionschemeisusedfor asfor all tracers; however, NICAM, ICON, …use different time stepsforcontinuityequationandtracers use time averagedmassfluxes! (talksbyH. Tomita, G. Zängl) M. Baldauf (DWD)
Semi-Lagrangianadvection (overviewgivenby M. Diamantakis) • asusedforthedynamicalcore in: • ARPEGE(MeteoFrance)/ • IFS(ECMWF)/ALADIN • UM(UKMO) • HIRLAM • SL-AV(Russia) • GEM(Environment Canada) • GFS(NCEP) • GSM(JMA) • … • Fortraceradvection, localconservationpropertiesare not easy toachieve: • e.g. SLICE-3D (Zerroukat, Allen (2012) QJRMS) implemented in ENDGameversionof UM (but computationally expensive) • on unstructuredgrids: conservativeSLprobablytoo expensive. • (conservative) SL schemes on structuredgridsmaybebeneficialformanytracers M. Baldauf (DWD)
Finite-element basedmethods • severalContinuous Galerkin (CG) / spectralelements /Discontinuous Galerkin (DG) developments: • NUMA (Navalpostgraduateschool, Monterey) F. Giraldo, … • HOMME (NCAR) R. Nair, … • CAM-SE (Sandia Nat. lab.) M. Taylor • COSMO (DWD) firststepsby D. Schuster, … • … • generalremarksin talk ofC. Cotter M. Baldauf (DWD)
C. Cotter (Imperial Coll.) linear 1D wave equation M. Baldauf (DWD)
C. Cotter (Imperial Coll.) Ladyzhenskaya/Babuska/Brezzi (LBB) condition guarantees compatibility M. Baldauf (DWD)
Discontinuous Galerkin (DG) methods in a nutshell talk by F. Giraldo • weak formulation(increases solution space) From Nair et al. (2011) in ‚Numerical techniques … Finite-element ingredient e.g. Cockburn, Shu (1989) Math. Comput. Cockburn et al. (1989) JCP e.g. Legendre-Polynomials Finite-volume ingredient Lax-Friedrichs flux ODE-systemforq(k)jl M. Baldauf (DWD)
1D wave expansion with a Discontinuous Galerkin (DG) discretization Literature: Hu, Hussaini, Rasetarinera (1999) JCP: 1D advection-, 2D wave-equation Hu, Atkins (2002) JCP: non-uniform grids k=k() Ainsworth (2004) JCP: multi-dim. advection equation M. Baldauf (DWD)
Numerical dispersion relation for the 1D wave equation Re x/c DG with p=0,1,2,3 (=c used) kx Im x/c • max || x/c 1 + 2.6 p + 0.33 p² • increases slightly stronger than linear with p. • Choose not too large p! kx M. Baldauf (DWD)
M. Baldauf (DWD) F. Giraldo (NPS)
Other remarks • „WGNE tableaboutcomputerressources 2014“uptonow, I havereceivedcontributionsfromChiashi (Japan), Hoon (Korea), Xueshun (China), Ayrton (Canada),Jean-Noel (ECMWF) and DWD • shallwecontinuethesurveyaboutdynamicalcores (by Mikhail Tolstyk) ? M. Baldauf (DWD)