A discontinuous Galerkin approach for atmospheric flows with implicit condensation

TL;DR

This work develops a high-order discontinuous Galerkin method for moist atmospheric flows with implicit condensation by aggregating water vapour and cloud water into a single moist density $\rho_m$ and recovering $\rho_v$ and $\rho_c$ via a pointwise nonlinear solve after each time step. The approach enables explicit time stepping with a diagonal mass matrix and a matrix-free implementation, enhanced by a local artificial diffusion stabilization and an explicit sponge layer to suppress top-boundary reflections. The method is validated through multiple 2D and 3D test cases, including inertia-gravity waves, Bryan–Fritsch moist benchmark, mountainous terrain, rising thermals, and a squall-line scenario, demonstrating optimal high-order convergence (up to $k+1$) in moist cases without rain and robust stability with rain. The results on structured and unstructured meshes, together with strong parallel scalability, indicate that the framework is well-suited for large-scale atmospheric simulations on HPC platforms, with open-source code provided for reproducibility.

Abstract

We present a discontinuous Galerkin method for moist atmospheric dynamics, with and without warm rain. By considering a combined density for water vapour and cloud water, we avoid the need to model and compute a source term for condensation. We recover the vapour and cloud densities by solving a pointwise non-linear problem each time step. Consequently, we enforce the requirement for the water vapour not to be supersaturated implicitly. Together with an explicit time-stepping scheme, the method is highly parallelisable and can utilise high-performance computing hardware. Furthermore, the discretisation works on structured and unstructured meshes in two and three spatial dimensions. We illustrate the performance of our approach using several test cases in two and three spatial dimensions. In the case of a smooth, exact solution, we illustrate the optimal higher-order convergence rates of the method.

Figures (14)

• Figure 1: Example 1: Convergence results of the primary and secondary variables, and wet equivalent potential temperature on a series of structured quadrilateral meshes for different polynomial orders.
• Figure 2: Example 1: Inertia gravity waves in a saturated atmosphere. Perturbation of the wet equivalent potential temperature from the hydrostatic base state. Left: Profile along the line $z=5\,\unit{km}$ at $t=2520\,\unit{s}$ for different polynomial orders on a structured quadrilateral mesh with $h=125\,\unit{m}$, centred around $x=200\,\unit{km}$. Right: Solution at $t=3600\,\unit{s}$ computed using $\mathbb{Q}^3$ elements with $h=500\,\unit{m}$, centred around $x=222\,\unit{km}$.
• Figure 3: Example 1: Strong scaling results for the time-stepping loop using NGSolve's shared memory (SM) and distributed memory (MPI) parallel implementations. The number of degrees of freedom (dofs) are counted for the finite element space for the eight primary variables.
• Figure 4: Example 2: Density potential temperature perturbation at $t=1000\,\unit{s}$ in the region $(6\,\unit{km}, 14\,\unit{km})\times(3\,\unit{km},8.5\,\unit{\km})$ on a series of different triangular and quadrilateral meshes using different order finite elements.
• Figure 5: Example 3: Cloud density ( ^2) and contours at $t=3600\,\unit{s}$ resulting from inertia gravity waves in a saturated atmosphere without initial clouds. Computed using $\mathbb{Q}^3$ elements on a structured quadrilateral mesh. Left: Profiles along the lines $z=5\,\unit{km}$ and $7.5\,\unit{km}$. Right: Density in the whole domain, computed with $h=500\,\unit{m},\Delta t=0.15\,\unit{s}$.

Theorems & Definitions (2)

• Remark 3.1
• Remark 3.2: Loss of higher-order convergence