630 likes | 802 Views
Michael Baldauf Deutscher Wetterdienst, Offenbach, Germany. Working group 2: Dynamics and Numerics report ‘Oct. 2006 – Sept. 2007’ COSMO General Meeting, Athen 18.-21.09.2007. 2.11 Alternative discretizations (due to alternative grids)
E N D
Michael Baldauf Deutscher Wetterdienst, Offenbach, Germany Working group 2:Dynamics and Numerics report ‘Oct. 2006 – Sept. 2007’COSMO General Meeting, Athen18.-21.09.2007 19.09.2007
2.11 Alternative discretizations (due to alternative grids) 2.11.1 [RENUMBERED] Remove grid redundancy by Serendipity GridsDWD: Steppeler 09/05 Report The serendipity grids should be investigated, which reduce the redundancy of the interpolation procedures. In this way they achieve more accuracy and more efficiency. 2.11.2 [NEW] Higher order discretization on unstructured grids using Discontinuous Galerkin methods DWD: Baldauf, Univ. Freiburg: Kroener, Dedner, NN. 2009 start, 2011 report In the DFG priority program 'METSTROEM' a new dynamical core for the COSMO-model will be developed. It will use Discontinuous Galerkin methods to achieve higher order, conservative discretizations. Currently the building of an adequate library is under development. The work with the COSMO-model will start probably at the end of 2009. This is therefore base research especially to clarify, if these methods can lead to efficient solvers for NWP. 19.09.2007
2.3.1 Radiative upper boundarycondition DWD: Herzog 09/05 Report The Klemp Durran boundary is further developed. 2.3.2 [NEW] Radiative upper boundary condition; non-local in time NN report in 06/2009 At the University Freiburg a Radiative upper boundary condition was developed. It is non-local in time, but nevertheless can be implemented efficiently into non-hydrostatic models. This radiation condition will be further developed during the DFG priority program METSTROEM. 19.09.2007
2.6.3 Implementation of neglected diabatic terms in p'-equation MPI-H: Petrik, DWD: Herzog 2.10 Diagnostic tools 2.10.1 [NEW] Application of the integration tool to energy, mass balance DWD: Baldauf, MPI-H: Petrik The integration tool to calculate balance equations by volume integrations of densities and surface integrations of fluxes developed in the Priority project 'Runge-Kutta', Task 3 will be applied to questions of energy and mass budgets. 19.09.2007
Michael Baldauf Deutscher Wetterdienst, Offenbach, Germany COSMO Priority Project:Further developments of the Runge-Kutta Time Integration Scheme report ‘Oct. 2006 – Sept. 2007’COSMO General Meeting, Athen18.-21.09.2007 19.09.2007
Tasks of the Priority Project ‚Runge-Kutta‘: • Looking at pressure bias • Continue RK case studies • Conservation • Advection of moisture quantities in conservation form • Investigation of convergence • Deep valleys • (Different filter options for orography) (finished) • Higher order discretization in the vertical for Runge Kutta scheme • Physics coupling scheme • Testing of alternative fast wave scheme • Development of a more conservative dynamics (planned) • Development of an efficient semi-implicit solver in combination with RK time integration scheme (planned) • NEW: Divergence damping in a truely 3D-version • NEW: DFI for RK 19.09.2007
List of people contributing to the project (Sept. 2006 - Aug. 2007): • (alphabetical order) • Michael Baldauf (DWD, D) • Gabriella Ceci (CIRA, I) • Guy deMorsier (MeteoCH, CH) • Jochen Förstner (DWD, D) finishes work in COSMO • Almut Gassmann(Univ. Bonn, D) finishes work in COSMO • Lucio Torrisi (CNMCA, I) • Pier Luigi Vitagliano (CIRA, I) • new: Gdaly Rivin (Roshydromet, RU) Additional meeting of PP-RK-Group during the LM-User-Workshop, Langen, 07./08.03.2007 19.09.2007
Task 1: Looking at pressure bias (Torrisi) Goals: verifications of LM 7 km runs showed a higher positive pressure bias for the RK core than for the Leapfrog core, whereas other variables show comparable behaviour. Reasons and solutions? RK starting point of the task: Leapfrog 19.09.2007
Task 1: • Status • the measurements • p'T'-dynamics • dynamical bottom boundary condition (DBBC) • brought a certain improvement • Intermediate question: • Why is RK sensitive to DBBC whereas Leapfrog is not? • LME parallel run at DWD shows even reduced pressure bias compared to Leapfrog in period '09.-25.02.07' (but is this a somewhat special period because there was mostly a high pressure situation in Europe ?) • actual status of the problem is not clear • no further work done • Plans: • Test idealised mountain flow with RK/Leapfrog --> information, if a pressure bias is caused by dynamics or by parameterizations (T. Davies proposal) (G. Ceci) • diabatic terms in pressure equation (R. Petrik (MPI-Hamburg), H. Herzog (retired) ) see WG2, Task 2.6.3 • determine mass conservation with conservation inspection tool from task 3 (Petrik, Baldauf)see WG2, Task 2.10.1 19.09.2007
Task 1: RMSE 19.09.2007
Task 1: BIAS 19.09.2007
Task 1: RMSE 19.09.2007
Task 1: BIAS 19.09.2007
k+1 -k-1 zk+1 - zk-1 k+1 -k-1 2( zk+1/2 - zk-1/2 ) Intermediate question: Why is RK sensitive to DBBC whereas Leapfrog is not? Possible answer: Discretization of vertical gradients in p‘-equation Leapfrog: RK: Both versions are consistent and even of the same order of accuracy. But the first version seems to better coincide with DBBC. With the LF approach a slight improvement in MSLP bias is found 19.09.2007
Task 2: Continue RK case studies (Torrisi, deMorsier) extensive verification of the other tasks • Status: • CNMCA: • COSMO-ITA tests with LF / RK / RK+Semi-Lagr. (results see task 4) • COSMO-ITA tests (RK+SL) with nudging and 3D-VAR interpolated initial state Work to do: verifications should be continued 19.09.2007
Status: MeteoSwiss: since ~15.08.2006 runs a 2.2 km version with assimilation cycle. Several unstable cases found in previous winter period (e.g. ‚13. Jan. 2004‘) most of them could be simulated with Semi-Lagrange Adv. for moisture variables Winter storms Kyrill ('18.01.2007') and Lothar ('26.12.1999') simulated with MeteoCH new pre-operational model chain (2.2 km and 6.6 km): new configuration: 1.) WRF-like RK3 used (instead of TVD-RK3) (as found at DWD for the Kyrill case)2.) Semi-Lagrange-Adv. for moisture (instead of Bott-scheme)3.) new level distribution especially in boundary layer (cures problems with TKE scheme) 19.09.2007
New level distribution compared to DWD-models New level distribution L60.2 in pre-operational COSMO-S2 now similar to COSMO-DE / IFS in the boundary layer. This prevents instabilities (2*dz structures for TKE) 19.09.2007
Storm Lothar, forecast, 24h wind gusts 2.2km TVD-RK32.2km 'WRF'-RK3 nrdtau=5, Dt=10 sec. nrdtau=5, Dt=20sec not possible with 20 sec. 19.09.2007
Task 3: Conservation (Baldauf) Tool for inspection of conservation properties will be developed. balance equation for scalar : temporal change sources / sinks flux divergence integration area = arbitrarily chosen cuboid (in the transformed grid, i.e. terrain-following) • Status: available in LM 3.23: • Subr. init_integral_3D: define cuboid (in the transformed grid!), prepare domain decomp. • Function integral_3D_total: calc. volume integral Vijk Vijk • Subr. surface_integral_total: calc. surface integrals V jijk * Aijk • prelimineary idealised tests were carried out • report finished; will be published in the next COSMO-Newsletter Nr. 7 (2007) Task is finished (Study of conservation properties will be continued in collaboration with MPI-Hamburg, see WG2 Task 2.10.1) 19.09.2007
Task 3: Weisman-Klemp (1982)-test case without physical parameterisation (only advection & condensation/evaporation) Semi-Lagrange-Adv. for qx with multiplicative filling x := (qv + qc ) total moisture mass M = x dV (Mn-Mn-1) / t Res. violation in moisture conservation (?) total surface flux timestep 19.09.2007
Task 3: Weisman-Klemp (1982)-test case with warmer bubble (10 K) without physical parameterisation, without Condensation/Evap. Semi-Lagrange-Adv. for qx with multiplicative filling x := (qv + qc ) Residuum 0 advection seems to be ‚conservative enough‘ total moisture mass M = x dV (Mn-Mn-1) / t possible reasons for conservation violation: saturation adjustment conserves specific mass (and specific energy) but not mass (and energy) itself ! Res. total surface flux timestep 19.09.2007
Task 4: Advection of moisture quantities in conservation form (Förstner, Baldauf) implementation of the Bott (1989)-scheme into the Courant-number independent advection algorithm for the moisture densities (Easter, 1993, Skamarock, 2004, 2006) Status: Task was finished in Sept. 2006 because implemented schemes (Bott-2, Bott-4) behaved well But in the meanwhile: stability problems occured in some cases revival of the task necessary! 19.09.2007
Problems found with Bott (1989)-scheme in the meanwhile: 1.) Directional splitting of the scheme: Parallel Marchuk-splitting of conservation equation for density can lead to a complete evacuation of cellsSolution: Easter (1993), Skamarock (2004, 2006), mass-consistent splitting 2.) Strang-splitting ( 'x-y-z' and 'z-y-x' in 2 time steps) produces 2*dt oscillations Solution: proper Strang-Splitting ('x-y-2z-y-x') in every time step solves the problem, but nearly doubles the computation time 3.) metric terms prevent the scheme to be properly mass conserving Schär–test case of an unconfined jet and ‚tracer=1‘ initialisation (remark: exact mass conservation is already violated by the 'flux-shifting' to make the Bott-scheme Courant-number independent) 19.09.2007
COSMO-ITA: RK+SL / RK+new Bott SL Bott RK+new Bott has a larger FBI for all precipitation thresholds than RK+SL(= COSMO-ITA operational run). Moreover, RK+new Bott has a deterioration in MSLP bias and RMSE after T+12h. 19.09.2007
Semi-Lagrangian advection in COSMO-model Moisture transport DWD: COSMO-DE: Bott-scheme used COSMO-EU: SL scheme planned operationally MeteoCH: COSMO-S2 and COSMO-S7: SL scheme used pre-operationally CNMCA: COSMO-ITA 2.8: SL-scheme used pre-operationally SL is not positive definite clipping necessary 'multiplicative filling' combines clipping with global conservation Work done: 'multiplicative filling uses global summation (summation of REAL not associative and therefore not reproducible) New subroutine sum_DDI( f(:,:) )(with accuracy estimation) reproducibility for 'multiplicative filling' is now reached 19.09.2007
Task 5: investigation of convergence • (Ceci, Vitagliano) • Goals: determination of the spatial and temporal order of convergence of the RK-scheme in combination with advection schemes of higher order. • Planned test cases: • linear, 2D, hydrostatic mountain flow (h=10 m, a=10 km) • linear, 2D, non-hydrostatic mountain flow (h=10 m, a=500 m) • nonlinear, 2D mountain flows (dry case) (h=500 m, a=10 km) • linear, 3D mountain flow • nonlinear mountain flows with precipitation 19.09.2007
Task 5: investigation of convergence • Work done: • use of LM 3.21 ( p'T'-dynamics) • Problems with unrealistic vertical winds at lateral boundaries solved:Introduction of an appropriate reference state p', T'=0 at the beginning • upper Rayleigh damping layer:idealised test cases need a properly chosen damping constant.The range for it reduces by decreasing the number of levels in damping layer( compare new WG2 task: 2.3.2, radiative upper boundary condition, non-local in time) • Calculation of L1, L – error norms for the dry 2D test cases • Additionally calculation of mountain drag coefficient, and vertical momentum flux • analytical solution for linear mountain flow (Klemp, Lilly, 1978) programmed • Plans: • determine L2, L – errors of KE, w, ..., dependent from x, t, ... for the test cases 19.09.2007
Task 5: Problems with unrealistic vertical winds at lateral boundaries solved Spurious vertical velocity at lateral boundaries original standard reference atmosphere and p*, T* modified to get constant Brunt-Väisälä frequency. comparison between the original and modified standard reference atmosphere. The effect gets larger as the time step decreases 19.09.2007
Task 5: investigation of convergence Comparison with analytical solution solution with a damping layer of 85 levels and nRΔt=200. Analytical solution following Klemp-Lilly (J.Atmos.Sc. 35, 78-107, 1978) 19.09.2007
CONVERGENCE OF VERTICAL VELOCITY w L1 = average of errors L = maximum error Convergence slightly less than 2. order. (2. order at smaller scales?) 19.09.2007
NON LINEAR HYDROSTATIC FLOW: STREAMLINES Stable and stationary solution of this non-linear case! 19.09.2007
CONVERGENCE OF VERTICAL VELOCITY w L1 = average of errors L = maximum error 19.09.2007
Task 6: deep valleys (Förstner, Torrisi, Reinhardt, deMorsier) Goal:detection of the reason for the unrealistic ‚cold pools‘ in Alpine valleys Status (old):The reason for the cold pools was identified: metric terms of the pressure gradient Dynamical Bottom boundary condition (DBBC) (A. Gassmann (2004), COSMO-Newsl.) and a slope-dependent orography-filtering cures the problem to a certain extent. • Proposal for future work: • inspect the limitations of the terrain following coordinate for steeper slopes, e.g. • for future LMK ~1 km horizontal resolution • for application of aLMo 2 (MeteoCH) in Alpine region • Does the strong conservation form (Gassmann-scheme) delivers advantages for steep terrain? 19.09.2007
Task 7: Different filter options for orography (Förstner) Status: the orography filtering is now sufficiently weak for DWD-LMK applications (max. slopes 30% allowed) finished 19.09.2007
Task 8: Higher order discretization in the vertical for RK-scheme (Baldauf) Improved vertical advection for the dynamic var. u, v, w, T (or T‘), p‘ motivation: resolved convection • vertical advection has increased importance => use scheme of higher order (compare: horizontal adv. from 2. order to 5. order) • => bigger w (~20 m/s) => Courant-crit. is violated =>implicit scheme or CNI-explicit scheme up to now: implicit (Crank-Nicholson) advection 2. order (centered differences) new: implicit (Crank-N.) advektion 3. order LES with 5-banddiagonal-matrix but: implicit adv. 3. order in every RK-substep; needs ~ 30% of total computational time! planned: use outside of RK-scheme (splitting-error?, stability with fast waves?) Status: no further work done Work to do: best combination with time integration scheme? 19.09.2007
Task 9: Physics coupling scheme (deMorsier, Förstner, Baldauf, NN) original idea: problems with reduced precipitation could be due to a nonadequate coupling between physics scheme and dynamics Problems in new physics-dynamics coupling (NPDC): • Negative feedback between NPDC and operational moistturbulence parameterization (not present in dry turbulence parameterization) • 2-z - structures in the specific cloud water field (qc) • 2-z - structures in the TKE field, unrealistic high values, where qc > 0 Status: no further work done • Work to do: • what are the reasons for the failure of the WRF-PD-scheme in LM? (turbulence scheme?) • Test different sequences of dynamics and physics (especially physics after dynamics) test tool (Bryan-Fritsch-case) is developed in PP ‚QPF‘, task 4.1 19.09.2007
Task 10: Testing of alternative fast wave scheme • (Gassmann, Torrisi, Förstner) • Goals: • p‘T‘-RK-scheme • ‚shortened-RK2‘-scheme (Gassmann) • this allows the use of the ‚radiative upper boundary condition‘ (RUBC) • Properties of A. Gassmann dyn. core: • Splitting up of vertical advection of p*/T into fast/slow mode equations and consistent boundary conditions • Vertical average to half levels: mass weighted mean (in RK simple mean) and base-state consistent formulation of the discrete w-equation • Different horizontal pressure gradient discretization • Divergence in conservative flux-form • Slightly different buoyancy term • No divergence damping 19.09.2007
Status: • The fast waves part (Gassmann) is combined with the Leapfrog scheme in LM 3.21 • Original Gassmann dynamical core poses stability problems in several cases! • Gassmann fast waves part in RK3 worked in only 1 case • ‚shortened RK2‘-scheme (Gassmann (2002), Gassmann and Herzog (2007)) is implemented into LM 3.21 using the fast waves solver of RK3 and the RK3 advection/physics subroutines • Preliminary investigation of this dynamical core (L. Torrisi)tested in real cases for a five days period: similar results to the RK3 splitting method • Separate inspection of divergence in conservation form and vertical staggering • Implementation questions pointed out: • Splitting of contravariant vertical velocityposes problems in formulation of lower boundary conditions • Work to do: • Further tests within CLM-community (dx=18 km) (B. Früh, Univ. Karlsruhe) 19.09.2007
Gassmann Time-Splitting Method 19.09.2007
Task 13: Divergence damping in a truely 3D-version (NEW) (Baldauf) Description: Cases occured, where the up to now used 'quasi-2D' divergence filtering lead to unstable results. But a complete abandoning of the divergence filtering (as proposed by A. Gassmann for her dynamical core) also leads to several instabilities. This was also shown by stability analyses of the RK-core by M. Baldauf. P. Prohl (DWD) could demonstrate, that the Bryan-Fritsch- test case of a rising warm bubble is unstable with 'quasi-2D' divergence damping but becomes stable only with a full 3D (=isotropic) version (realised with a preliminary explicit formulation). For operational use an implicit version of 3D divergence damping is necessary. (this task was defined as a subtask in Task 10; but Task 10 concentrates more and more to the A. Gassmann dynamical core) 19.09.2007
Task 14: DFI for RK NEW (L. Torrisi) Description: Bug of DFI has been repaired in the Leapfrog-version of the COSMO-model. This has still to be done in RK (L. Torrisi already investigated some work) 19.09.2007
List of people contributing to the WG2 (‚Numerics‘): • (alphabetical order) • Euripides Avgoustoglou (GR) LM-Z • Michael Baldauf (DWD, D) • Heinz-Werner Bitzer (MetBW, D) LM-Z, but mostly Meteorol. Analysis • Gabriella Ceci (CIRA, I) • Guy deMorsier (MeteoCH, CH) • Jochen Förstner (DWD, D) (has left LMK-comm. ICON) • Almut Gassmann (MPI, Hamburg) (has left LMK-comm. ICON) • Hans Herzog retired • Gdaly Rivin (Roshyd., RU) • Uli Schättler (DWD, D) • Jürgen Steppeler retired • Lucio Torrisi (CNMCA, I) • Pier Luigi Vitagliano (CIRA, I) ‚Numerics‘-group has reduced personal in future! 19.09.2007
Task 1: • Status (until Sept. 2006): • 5-day verifications were done for several model configurations • only little impact on PMSL by: • physics coupling • advection of qx (Bott, Semi-Lagrange) • most significant impact on PMSL: • new dynamical bottom boundary condition (DBBC) ( Task 6) • p‘T‘-dynamics in the RK-core • both measurements reduce the pressure bias verification area (L. Torrisi) 19.09.2007
Task 1: positive impact of p‘T‘-dynamics 19.09.2007
Task 1: positive impact of DBBC 19.09.2007
Task 1: Status: verifications carried out Work to do: pressure bias improved by another fast waves solver?? ( --> task 10) 19.09.2007
PP RK 2.1.4 + 2.1.3 Transport of Tracerin a Real Case Flow Field Bott (2nd)“Flux Form- DIV”+ Clipping init semi-Lagrange(tri-cubic)+ Clipping Bott (2nd)“Conserv. Form” 19.09.2007