Efficient and Conservative Approaches for COSMO-GM Model

1. Priority Project CDC Task 2: The compressible approach COSMO-GM, 06.-10.09.2010, Moscow Pier Luigi Vitagliano (CIRA), Michael Baldauf (DWD)

2. Task 2.3: Fully 3D, i.e. non-direction splitted, conservative advection scheme At MeteoCH, a diploma thesis (M. Müllner) has been started, toimplement MPDATA into COSMO. No serious problems expected.

3. Task 2.2: Complete Finite-Volume solver for the EULER equations MOTIVATIONS AND GOALS • Improve numerical efficiency • Improve conservation properties • Capability to deal with steep orography • Test a new time integration scheme • Test spatial schemes based on finite volumes

4. Dry Euler equations (without Coriolis force) in conservative form W= Fx = Fy = Fz = B= E = Ekin + Eint = ½ U2 + cp T

5. SPATIAL DISCRETISATION • Finite Volumes approach • Integral form allows discontinuities in the flow field • Conservation laws applied to each sub-domain (cell) • Variables stored at cell centers • Fluxes approximated at cell face centers (VW)/t + R(W) = 0 R(W) = Q – B – D Q = fluxes B = source terms D = k∆4artificial dissipation

6. DUAL TIME STEPPING Wn+1/t + ½(3Wn+1- 4Wn + Wn-1)/Dt + R(Wn+1) = 0 formulation is A-stable and damps the highest frequency  very large physical time step Dt can be used • Solution of the implicit equation system: • add a pseudo-time tderivative to the unsteadyequation • integration in  is performed by an explicit Runge-Kutte scheme • advance the solution in t until the residual of the unsteady equation is negligible • convergence acceleration techniques can be adopted without loss of time accuracy:residual averaging, local time stepping, multigrid Jameson, A., 1991: Time Dependent Calculations Using Multigrid,with Applications to Unsteady Flows Past Airfoils and Wings. AIAA Paper 91–1596

7. DUAL TIME STEPPING Example of time integration with DTS: a norm of the residuals of mass transport equations is monitored

8. P·W/t + R(W) = 0 PRECONDITIONING Improve convergency in dual time for low Mach number flows Correct ill-behaved artificial viscosity fluxes at low Mach Difficulties rise from large ratio between acoustic wave speed and fluid speed Premultiplying the time derivative changes the eigenvalues of the system and accelerates the convergence to steady state. Turkel, E., 1999: Preconditioning techniques in computational fluid dynamics. Annu.Rev.Fluid Mech. 1999,31:385-416. Venkateswaran, S., P. E. O. Buelow, C. L. Merkle, 1997: Development of linearized preconditioning methods for enhancing robustness and efficiency of Euler and Navier-Stokes Computations, AIAA Paper 97-2030.

9. PRECONDITIONING Example of convergence to steady solution with and without Preconditioning

10. For task 2.2 the following idealised dry test cases where defined: (only a reduced set of the test cases from task 3.1) • Atmosphere at rest (G. Zaengl (2004) MetZ) • test balance of pressure gradient forces, metric correction terms and buoyancy • Cold bubble (Straka et al. (1993)) strong nonlinear, unstationary test, well established reference solution • Mountain flow tests:stationary test, well known (partly analytic) solutions • Schaer et al (2002) sect. 5b • Bonaventura (2000) JCP • Linear Gravity waves (Skamarock, Klemp (1994), Giraldo (2008))unstationary test, wave expansion, analytic solution available

11. ATMOSPHERE AT REST Pressure gradient discretisation Field initialisation Effect of mesh skewness Flux – force unbalance standard formulation:pmk= ½ (pm+ pk) pressure gradient correction: pmk= ½ {pm + (∂p/∂z)mΔzm+ pk + (∂p/∂z)k Δzk} (∂p/∂z)m= ρmg

12. ATMOSPHERE AT REST INITIAL FIELD W component

13. ATMOSPHERE AT REST standard pressure gradient pressure gradient correction W component after 90000 seconds Δt=3 sec

14. ATMOSPHERE AT REST Solution with pressure gradient correction after 90000 seconds Δt=300 sec

15. ATMOSPHERE AT REST CONCLUSIONS Solution is not affected by physical time step nor by CFL Pressure gradient correction improves initialisation of the fields, but has only a smaller positive influence on accuracy after 1 day simulation (but helps convergence).The induced vertical velocities are in the same order of magnitude than in COSMO. Some issues with boundary conditions

16. TEST CASE MOUNTAIN FLOW Linear, hydrostatic case Flow over a gaussian mountain simulated with a test code based on finite volumes conservative schemes. Vertical velocity component. The dashed line shows the lower boundary of the Rayleigh damping layer, which prevents the wave reflection.

17. STRAKA TEST • Implemented viscous fluxes with constant ν • Implemented reference atmosphere with constant ∂T/∂z INITIAL FIELD

18. STRAKA TEST Δt=10 sec SOLUTION AFTER 600 sec Δx=50 m reference solution by Straka et al (1993) 4.8 km 19.2 km

19. STRAKA TEST • good agreement with reference solution by Straka et al. (1993) • time step from 0.25 to 100 sec possible • mesh size from 25 m to 200 m • solution diffused with larger time steps

20. Gravity wave test (Skamarock, Klemp (1994) MWR)

21. Conclusions • most of the test cases were carried out successfully • atm. at rest • cold bubble • linear, hydrostatic mountain flow • linear gravity wave test • some idealised tests are still missing • not all mountain flow tests availablemost probably due to initialisation/setup problems of the test cases • a lot of work had to be done with the implementation of buoyancy termsinto the model; this looks promising • dual time stepping is a promising time integration approach;statements about efficiency compared to current COSMO are not so easy until now; but scalability is probably not an obstacle

22. Scalability on future supercomputing platforms • no tests made yet (toy model) • a 3D implicit solver is used, but in an 'explicit' manner due to the dual time stepping  should not pose serious problems • preconditioning: only local operations  linear speedup expected(Choi, Merkle (1993) JCP)

23. Task 2.2.2: It has to be clarified how the moist equations should be formulated. An adequate test case should be performed: Weisman, Klemp (1982) (warm , moist bubble test) Task 2.2.3: The properties of the A-grid formulation concerning wave propagation should be investigated. Therefore a wave analysis on the A-grid will be performed.

24. Task 2.5: Testing the dual time stepping in COSMO The dual time stepping (DTS) scheme is generally able to integrate implicit equation systems. Therefore it can be used to integrate the COSMO equations by abandoning the time-splitting procedure. Fast processes have to be formulated implicitely, but with the same spatial discretizations as they are used now for the Runge-Kutta scheme. Such an implementation is not expected to be more efficient, but possibly could solve problems connected with steep terrain. This task therefore serves as an intermediate step towards task 2.6 This preliminary DTS implementation can be at first tested with the implemented idealised test cases (see task 3.1). This testing can be performed ‘by a press of a button’ in COSMO. Deliverables: COSMO model using DTS scheme.

25. Task 2.6: Implement the Finite Volume solver into COSMO Finally the scheme developed in task 2.2 will be implemented into COSMO. Again with the implemented idealised test cases (see task 3.1) a testing to find elementary bugs can be performed ‘by a press of a button’ in COSMO. Deliverables: COSMO model using the compressible, implicit Finite Volume dyn. core Task 2.7: Perform realistic test cases After finishing task 2.6 real case simulations with full physics parameterisations with COSMO are possible. Stand-alone runs for several weather regimes can be performed for both dynamical cores (FV, RK) at different resolutions. One has to obey that physical parameterizations have to be adapted to the new dynamical core. This probably requires support from the physical parameterization working group. Deliverables: report about the behaviour of real case test simulations with COSMO