A semi-analytical pseudo-spectral method for 3D Boussinesq equations of rotating, stratified flows in unbounded cylindrical domains

TL;DR

This paper introduces a high-accuracy spectral method with an exponential time integrator for simulating 3D rotating, stratified flows in unbounded cylindrical domains, effectively handling fast wave dynamics and background shear.

Contribution

It develops a novel semi-analytical pseudo-spectral approach combined with an ETD scheme that overcomes stability constraints in simulating complex geophysical and astrophysical flows.

Findings

Accurate conservation of energy and angular momentum demonstrated.

Significant performance improvements over traditional schemes achieved.

Validated stability and accuracy for complex flow simulations.

Abstract

We present a pseudo-spectral method for solving the three-dimensional Boussinesq equations in unbounded cylindrical domains, specifically tailored for rotating, stably stratified flows subject to strong azimuthal shear. To effectively capture the global geometry without sacrificing spectral accuracy, the spatial discretization employs Fourier expansions in the azimuthal and axial directions alongside mapped associated Legendre polynomials in the radial direction. This basis spans the semi-infinite domain while analytically resolving the coordinate singularity at the origin. While this spectral framework ensures high spatial fidelity, the temporal integration of these rotating shear flows presents a formidable computational challenge due to the numerical stiffness driven by fast restorative wave forces and rapid background advection. To circumvent this, we develop an exponential time differencing (ETD) scheme that analytically integrates the fully coupled linear operator, including the radially dependent advective cross terms. By encoding the physical resonance characteristics and stability limits of the background flow directly into the integration operators, the proposed ETD formulation removes the numerical stability constraints imposed by the background shear and stratification. This permits integration time steps scaled by the slow macroscopic evolution of the physical instabilities rather than the fast background kinematics, offering significant performance gains over standard mixed implicit-explicit schemes. The method's accuracy and stability are validated through the precise conservation of energy and angular momentum, establishing a robust framework for simulating instabilities in astrophysical and geophysical vortices.