Spurious Diapycnal Diffusion in Terrain-following Coordinates : New Challenges

1. Noumea, New Caledonia Spurious Diapycnal Diffusion in Terrain-following Coordinates : New Challenges P. Marchesiello Contributions from X. Couvelard, L. Debreu & C. Menkes 2007 ROMS/TOMS Workshop Meeting

2. Context: Workshop on Numerical Methods in Ocean Models, Bergen, August 2007 • G. Danabasoglu: “The sigma-coordinate models are not at a disadvantage for climate applications” • Pressure Gradient Errors: “New methods, i.e., better numerics, reduce the PGF error considerably” • r ~ 0.5 !!! • Spurious Diapycnal Mixing: « Because of better advection schemes and higher resolution, spurious diapycnal mixing errors are not significant in regional models. These errors are radiated away through open boundaries. » • S. Griffies: “Previous generation of advection schemes [Quick, FCT] commonly used in [geopotential] ocean models exhibit unacceptable levels of diapycnal mixing” • K~1-10 cm²/s in eddy-resolving simulations Are we overly optimistic? 2007 ROMS/TOMS Workshop Meeting

3. Spurious Mixing: Evidences Salinity at 1000m, characteristics of AAIW 2007 ROMS/TOMS Workshop Meeting

4. Tracer transport in the ocean • In a stratified ocean, turbulent tracer transport is strongly anisotropic at the scale above 3D turbulent scale (10-100m): AD=0.1 cm²s-1 AI=107 cm²s-1 • Tracer transport in ocean models appears as: • Resolved advection processes • Parameterized subgrid scale processes • Physical closure: unresolved eddy advection. Use of isopycnal explicit mixing operators. • Numerical closure: clean up noise, avoiding unphysical mixing 2007 ROMS/TOMS Workshop Meeting

5. Numerical closure:Basic review of diffusive advection schemes • The tracer advection equation should be non-dispersive: the phase speed ω/k and group speed δω/δk are equal (=U) • The discretized equation using second order centered scheme is dispersive. We need to solve an advection-diffusion equation: • Condition of stability (Bryan et al., 1975) • Pe < 2, where Pe = UΔx/A • A = ½UΔx : equivalent to a first order upstream advection scheme which is monotonic (good) but too diffusive • Solutions (Shchepetkin & McWilliams, 1998): • Hyperdiffusion for better scale selection • But what value for B? • Diffusive high-order advection schemes • Oscillatory schemes (linear): Quick, UP3 • Monotone schemes (nonlinear) using antidiffusion with limiters: FCT, TVD, Super-B, Sweby, MDPPM, Prather … 2007 ROMS/TOMS Workshop Meeting

6. Upstream 1rst & 3rd order Centered 2nd & 4th order 2nd order with «limiters » 3rd order with «limiters » Durran, 2004 2007 ROMS/TOMS Workshop Meeting

7. The UP3 advection scheme and its split version • Quick-type (Leonard, 1979) or UP3 scheme: • Quadratic upstream interpolation of the tracer field to cell interfaces • 1/6 (rather than 1/8) removes O(δx²) error • Split-UP3: Holland et al. (1998) for NCOM and Webb et al. (1998) for OCCAM (their goal was to ensure stability using leapfrog): • Third-order upstream scheme is equivalent to: • 4th-order centered advection scheme • Biharmonique diffusion • Variable B Cubic polynomial interpolation 2007 ROMS/TOMS Workshop Meeting

8. Scaling analysis of spurious mixing using UP3 in sigma models • Redi (1982): • The maximum slope S is at the bottom: • Lee et al. (2002): • if tracer variance is associated with a stream of eddies, characterized by Rossby radius L, we can write: • Rossby radius L0 • If Δx < L0 , L = L0 • If Δx > L0 , L = Δx • Vertical diffusion scales as: Coarse High Coarse Medium High 2007 ROMS/TOMS Workshop Meeting

9. AD 1000m ? ΔS1000m Obtained from SWP ROMS compared with Levitus Scaling analysis and comparison with ROMS salinity errors Sraw = 10 % r0 = 0.2 U = 10 cm/s 2007 ROMS/TOMS Workshop Meeting

10. Solution: Split-UP3 with rotated mixing tensor • Rotation of biharmonic operator • Interpolation of diffusivity on cell interfaces • Need clipping of B because of vertical fluxes which are conditionally stable only • Problem is worst for isopycnal diffusion: need implicit soving? 2007 ROMS/TOMS Workshop Meeting

11. ROMS in the southwest Pacific 2007 ROMS/TOMS Workshop Meeting

12. Satellite / In-situ data Global ocean model Global atmospheric model WRF OpenDaP ROMS_TOOLS 3DVAR Initial & boundary conditions ROMS 2007 ROMS/TOMS Workshop Meeting

13. Satellite / In-situ data Global ocean model Global atmospheric model WRF OpenDaP ROMS_TOOLS 3DVAR Initial & boundary conditions ROMS 2007 ROMS/TOMS Workshop Meeting

14. Configuration of the southwest Pacific • ROMS Configuration (Marchesiello et al. 2003): • Climatological forcings • Horizontal resolution 1/12°, 1/6°, ¼°, ½° • Vertical resolution: 35 levels with θS=6, θB=0 • Run on a PC cluster (20 processors) 2007 ROMS/TOMS Workshop Meeting

15. Testing in a smaller configuration at ¼° resolution Topographic Slope S 2007 ROMS/TOMS Workshop Meeting

16. Split-UP3 along sigma surfaces LEVITUS SUP3 S The separation of terms in UP3 gives results that are similar to the original scheme, although the splitted version shows slightly larger diffusion UP3 2007 ROMS/TOMS Workshop Meeting

17. Split-UP3 along geopotential surfaces LEVITUS SUP3 The splitted version of UP3 is able to properly correct the diapycnal mixing problem UP3 2007 ROMS/TOMS Workshop Meeting

18. Split-UP3 with diffusion of tracer perturbations ( Mellor and Blumberg, 1985) Diffusion of tracer perturbations (instead of tracer fields) is performed along sigma surfaces (Mellor & Blumberg, 1985): The method gives reasonable results but introduces noise 2007 ROMS/TOMS Workshop Meeting

19. Hyperdiffusion Coefficient 1011 m²/s 2007 ROMS/TOMS Workshop Meeting

20. Southwest Pacific Aplication:Requirement of performances • Maintain the right water masses • Allow realistic eddy dynamics and energy spectrum (weak dissipation of potential energy) • Avoid numerical instability (clean up noise by respecting Peclet constraint) 2007 ROMS/TOMS Workshop Meeting

21. Maintaining water masses 2007 ROMS/TOMS Workshop Meeting

22. Allowing full energy spectrum (1/6°) UP3 Harmonique Diffusion SUP3 2007 ROMS/TOMS Workshop Meeting

23. Avoiding numerical instability Laplacian diffusion with Smagorinsky coefficient (4th order centered advection) SUP3 2007 ROMS/TOMS Workshop Meeting

24. Revisiting Smagorinsky and Peclet constraint • Smagorinsky (1963) coefficient for Laplacian diffusion • Griffies et al. (2000): If the deformation rate scales as: • CA ~ 0.1 (physical arguments): does not satisfy the Peclet constraint: CA > 0.5 ! This noise may be suppressed by excessive temporal filtering (Mellor et al., 1998) • Griffies et al. (2000): Smagorinsky-like coefficient for biharmonique diffusion • Theoretical value: CB=1/16 (from Peclet constraint) • Empirical value given: CB = 0.1-0.2 ~ twice the theoretical value D is deformation rate DT: tension and DS: shearing strain then C=0.5 C=0.2 Durran, 1991 2007 ROMS/TOMS Workshop Meeting

25. Conclusion • Large spurious diapycnal mixing arise in sigma models when using diffusive advection schemes • even the newest schemes are not likely to change that • not all implicit schemes can be split. • Can we devise advection schemes with built-in rotated diffusion? • Should we use temporal damping (not involving spurious spatial mixing) to overcome the Peclet constraint? • Split-UP3 is a solution that appears satisfying, but there is a bit of work to do: • Equivalence between original and split versions? • Clipping of diffusion because of numerical stability constraint • Implicit solving? • Isopycnal diffusion will be needed for large-scale problems, which complicates further the clipping problem • Need sorting out clipping problem 2007 ROMS/TOMS Workshop Meeting