Using nonlinear interactions to quantify and improve timestepping accuracy in the Rotating Shallow Water Equations

Abstract: A key objective for the next generation of weather and climate models is the reduction of wall-clock time. This is to ensure that predictions are made within operational timeframes and that power expenditure is minimised. One way to reduce the wall-clock time is by using a larger timestep; however, this can reduce the accuracy of a numerical solution. The extent of this reduction in solution quality can be difficult to ascertain, given the complexity of the oscillatory, multiscale, and nonlinear governing equations for geophysical fluid flow. This work details new investigations to quantify and improve the accuracy of large timestep solutions to geophysical fluid systems. We do so by focusing on the Rotating Shallow Water Equations (RSWEs). These are a simplified version of the full governing equations, yet retain the key numerical challenges of fast oscillations and multiple timescales. Our new research with the RSWEs can be categorised into two main topics: I: Quantification of nonlinear timestepping error through a new metric and test cases II: Improving large timestep accuracy of the phase-averaging method through a mean correction The novelty in this work is that we will investigate both topics from the perspective of nonlinear interactions and their role in the time evolution of the RSWEs. This allows us to examine how well nonlinear dynamics, which are constructed by these interactions, are replicated in a numerical model when taking a large timestep. For topic I of quantifying timestepping error, we reformulate the f-plane and doubly periodic RSWEs into a time evolution equation for spectral coefficients. This represents the dynamics as the result of nonlinear interactions of three linear waves in a triad. Linear wave dispersion errors can appear in the numerical representation of the triadic interaction frequencies and thus impact the nonlinear dynamics in a timestepped solution. We quantify this through a new triadic error, which combines stability polynomials from the oscillatory Dahlquist test equation. Advantageously, the triadic error reflects how discrepancies in modelling fast linear waves can impact the representation of slow nonlinear dynamics. Two new test cases for the RSWEs, including a class of triadic tests, are also developed for the analysis of time discretisations in complex numerical models. These cases focus on how accurately a large timestep solution replicates the slowly developing nonlinear dynamics amidst the fast oscillations. We demonstrate these test cases in three numerical models, including Gusto, and LFRic from the Met Office. Topic II introduces a new mean corrected version of the phase-averaging timestepping technique. Phase-averaging exploits the oscillatory nature of geophysical fluid equations to take larger timesteps, by mapping to a modulation variable and averaging the resulting timestepping function. The mean corrected version uses a modified mapping that includes averaged information about the nonlinear oscillations. This allows for reductions in averaging error and greater large timestep accuracy. The mean correction was introduced by Tao in [Commun Nonlinear Sci Numer Simulat 71 (2019) 1–21] for use in weakly nonlinear systems and with classical time averaging (over an infinite time period). We reformulate this technique for i) systems where the nonlinearity is not weak but the linear oscillations are fast, and ii) finite averaging windows instead of classical time averages. In particular, we introduce a local mean correction that combines Tao’s mean correction with the use of finite averaging; this retains low-frequency nonlinear oscillations that are removed in the classical mean correction. The mean corrected method is shown to improve phase-averaging accuracy in three numerical experiments, including the one-dimensional RSWEs for different linear oscillation speeds. In this test, we observe a best finite averaging window for computing the local mean correction. Finally, to investigate the potential of the mean correction in future weather and climate models, we integrate our new timestepping method into the phase-averaged finite element model of Yamazaki et al. [Q. J. R. Meteorol. Soc., 149 (755), 2504–2513, (2023)]. The code implementation required overcoming technical challenges, including the use of a Lawson timestepping algorithm and the computation of a regularised pseudoinverse. With resolutions to these, we present the first results of mean corrected phase-averaging for the RSWEs on the sphere.