Francis X. Giraldo

A study of spectral element and discontinuous Galerkin methods for the Navier-Stokes equations in nonhydrostatic mesoscale atmospheric modeling: Equation sets and test cases

We present spectral element (SE) and discontinuous Galerkin (DG) solutions of the Euler and compressible Navier-Stokes (NS) equations for stratified fluid flow which are of importance in nonhydrostatic mesoscale atmospheric modeling. We study three different forms of the governing equations using seven test cases. Three test cases involve flow over mountains which require the implementation of non-reflecting boundary conditions, while one test requires viscous terms (density current). Including viscous stresses into finite difference, finite element, or spectral element models poses no additional challenges; however, including these terms to either finite volume or discontinuous Galerkin models requires the introduction of additional machinery because these methods were originally designed for first-order operators. We use the local discontinuous Galerkin method to overcome this obstacle. The seven test cases show that all of our models yield good results. The main conclusion is that equation set 1 (non-conservation form) does not perform as well as sets 2 and 3 (conservation forms). For the density current (viscous), the SE and DG models using set 3 (mass and total energy) give less dissipative results than the other equation sets; based on these results we recommend set 3 for the development of future multiscale research codes. In addition, the fact that set 3 conserves both mass and energy up to machine precision motives us to pursue this equation set for the development of future mesoscale models. For the bubble and mountain tests, the DG models performed better. Based on these results and due to its conservation properties we recommend the DG method. In the worst case scenario, the DG models are 50% slower than the non-conservative SE models. In the best case scenario, the DG models are just as efficient as the conservative SE models.

A scalable spectral element eulerian atmospheric model (SEE-AM) for NWP: Dynamical core tests

Abstract A new dynamical core for numerical weather prediction (NWP) based on the spectral element method is presented. This paper represents a departure from previously published work on solving the atmospheric primitive equations in that the horizontal operators are all written, discretized, and solved in 3D Cartesian space. The advantages of using Cartesian space are that the pole singularity that plagues the equations in spherical coordinates disappears; any grid can be used, including latitude–longitude, icosahedral, hexahedral, and adaptive unstructured grids; and the conversion to a semi-Lagrangian formulation is easily achieved. The main advantage of using the spectral element method is that the horizontal operators can be approximated by local high-order elements while scaling efficiently on distributed-memory computers. In order to validate the 3D global atmospheric spectral element model, results are presented for seven test cases: three barotropic tests that confirm the exponential accuracy of...

Continuous and discontinuous Galerkin methods for a scalable three-dimensional nonhydrostatic atmospheric model: Limited-area mode

This paper describes a unified, element based Galerkin (EBG) framework for a three-dimensional, nonhydrostatic model for the atmosphere. In general, EBG methods possess high-order accuracy, geometric flexibility, excellent dispersion properties and good scalability. Our nonhydrostatic model, based on the compressible Euler equations, is appropriate for both limited-area and global atmospheric simulations. Both a continuous Galerkin (CG), or spectral element, and discontinuous Galerkin (DG) model are considered using hexahedral elements. The formulation is suitable for both global and limited-area atmospheric modeling, although we restrict our attention to 3D limited-area phenomena in this study; global atmospheric simulations will be presented in a follow-up paper. Domain decomposition and communication algorithms used by both our CG and DG models are presented. The communication volume and exchange algorithms for CG and DG are compared and contrasted. Numerical verification of the model was performed using two test cases: flow past a 3D mountain and buoyant convection of a bubble in a neutral atmosphere; these tests indicate that both CG and DG can simulate the necessary physics of dry atmospheric dynamics. Scalability of both methods is shown up to 8192 CPU cores, with near ideal scaling for DG up to 32,768 cores.

Implicit-Explicit Formulations of a Three-Dimensional Nonhydrostatic Unified Model of the Atmosphere (NUMA)

We derive an implicit-explicit (IMEX) formalism for the three-dimensional (3D) Euler equations that allow a unified representation of various nonhydrostatic flow regimes, including cloud resolving and mesoscale (flow in a 3D Cartesian domain) as well as global regimes (flow in spherical geometries). This general IMEX formalism admits numerous types of methods including single-stage multistep methods (e.g., Adams methods and backward difference formulas) and multistage single-step methods (e.g., additive Runge--Kutta methods). The significance of this result is that it allows a numerical model to reuse the same machinery for all classes of time-integration methods described in this work. We also derive two classes of IMEX methods, one-dimensional and 3D, and show that they achieve their expected theoretical rates of convergence regardless of the geometry (e.g., 3D box or sphere) and introduce a new second-order IMEX Runge--Kutta method that performs better than the other second-order methods considered. We...

A Conservative Discontinuous Galerkin Semi-Implicit Formulation for the Navier-Stokes Equations in Nonhydrostatic Mesoscale Modeling

A discontinuous Galerkin (DG) finite element formulation is proposed for the solution of the compressible Navier-Stokes equations for a vertically stratified fluid, which are of interest in mesoscale nonhydrostatic atmospheric modeling. The resulting scheme naturally ensures conservation of mass, momentum, and energy. A semi-implicit time-integration approach is adopted to improve the efficiency of the scheme with respect to the explicit Runge-Kutta time integration strategies usually employed in the context of DG formulations. A method is also presented to reformulate the resulting linear system as a pseudo-Helmholtz problem. In doing this, we obtain a DG discretization closely related to those proposed for the solution of elliptic problems, and we show how to take advantage of the numerical integration rules (required in all DG methods for the area and flux integrals) to increase the efficiency of the solution algorithm. The resulting numerical formulation is then validated on a collection of classical two-dimensional test cases, including density driven flows and mountain wave simulations. The performance analysis shows that the semi-implicit method is, indeed, superior to explicit methods and that the pseudo-Helmholtz formulation yields further efficiency improvements.

A spectral element shallow water model on spherical geodesic grids

SUMMARY The spectral element method for the two-dimensional shallow water equations on the sphere is presented. The equations are written in conservation form and the domains are discretized using quadrilateral elements obtained from the generalized icosahedral grid introduced previously (Giraldo FX. Lagrange‐ Galerkin methods on spherical geodesic grids: the shallow water equations. Journal of Computational Physics 2000; 160: 336‐368). The equations are written in Cartesian co-ordinates that introduce an additional momentum equation, but the pole singularities disappear. This paper represents a departure from previously published work on solving the shallow water equations on the sphere in that the equations are all written, discretized, and solved in three-dimensional Cartesian space. Because the equations are written in a three-dimensional Cartesian co-ordinate system, the algorithm simplifies into the integration of surface elements on the sphere from the fully three-dimensional equations. A mapping (Song Ch, Wolf JP. The scaled boundary finite element method—alias consistent infinitesimal finite element cell method—for diffusion. International Journal for Numerical Methods in Engineering 1999; 45: 1403‐1431) which simplifies these computations is described and is shown to contain the Eulerian version of the method introduced previously by Giraldo (Journal of Computational Physics 2000; 160: 336‐368) for the special case of triangular elements. The significance of this mapping is that although the equations are written in Cartesian co-ordinates, the mapping takes into account the curvature of the high-order spectral elements, thereby allowing the elements to lie entirely on the surface of the sphere. In addition, using this mapping simplifies all of the three-dimensional spectral-type finite element surface integrals because any of the typical two-dimensional planar finite element or spectral element basis functions found in any textbook (for example, Huebner et al. The Finite Element Method for Engineers. Wiley, New York, 1995; Karniadakis GE, Sherwin SJ. Spectral:hp Element Methods for CFD. Oxford University ( ˘

High-order triangle-based discontinuous Galerkin methods for hyperbolic equations on a rotating sphere

High-order triangle-based discontinuous Galerkin (DG) methods for hyperbolic equations on a rotating sphere are presented. The DG method can be characterized as the fusion of finite elements with finite volumes. This DG formulation uses high-order Lagrange polynomials on the triangle using nodal sets up to 15th order. The finite element-type area integrals are evaluated using order 2N Gauss cubature rules. This leads to a full mass matrix which, unlike for continuous Galerkin (CG) methods such as the spectral element (SE) method presented in Giraldo and Warburton [A nodal triangle-based spectral element method for the shallow water equations on the sphere, J. Comput. Phys. 207 (2005) 129-150], is small, local and efficient to invert. Two types of finite volume-type flux integrals are studied: a set based on Gauss-Lobatto quadrature points (order 2N-1) and a set based on Gauss quadrature points (order 2N). Furthermore, we explore conservation and advection forms as well as strong and weak forms. Seven test cases are used to compare the different methods including some with scale contractions and shock waves. All three strong forms performed extremely well with the strong conservation form with 2N integration being the most accurate of the four DG methods studied. The strong advection form with 2N integration performed extremely well even for flows with shock waves. The strong conservation form with 2N-1 integration yielded results almost as good as those with 2N while being less expensive. All the DG methods performed better than the SE method for almost all the test cases, especially for those with strong discontinuities. Finally, the DG methods required less computing time than the SE method due to the local nature of the mass matrix.

A spectral element semi-Lagrangian (SESL) method for the spherical shallow water equations

A spectral element semi-Lagrangian (SESL) method for the shallow water equations on the sphere is presented. The sphere is discretized using a hexahedral grid although any grid imaginable can be used as long as it is comprised of quadrilaterals. The equations are written in Cartesian coordinates to eliminate the pole singularity which plagues the equations in spherical coordinates. In a previous paper [Int. J. Numer. Methods Fluids 35 (2001) 869] we showed how to construct an explicit Eulerian spectral element (SE) model on the sphere; we now extend this work to a semi-Lagrangian formulation. The novelty of the Lagrangian formulation presented is that the high order SE basis functions are used as the interpolation functions for evaluating the values at the Lagrangian departure points. This makes the method not only high order accurate but quite general and thus applicable to unstructured grids and portable to distributed memory computers. The equations are discretized fully implicitly in time in order to avoid having to interpolate derivatives at departure points. By incorporating the Coriolis terms into the Lagrangian derivative, the block LU decomposition of the equations results in a symmetric positive-definite pseudo-Helmholtz operator which we solve using the generalized minimum residual method (GMRES) with a fast projection method [Comput. Methods Appl. Mech. Eng. 163 (1998) 193]. Results for eight test cases are presented to confirm the accuracy and stability of the method. These results show that SESL yields the same accuracy as an Eulerian spectral element semi-implicit (SESI) while allowing for time-steps 10 times as large and being up to 70% more efficient.

Semi-Implicit Formulations of the Navier-Stokes Equations: Application to Nonhydrostatic Atmospheric Modeling

We present semi-implicit (implicit-explicit) formulations of the compressible Navier-Stokes equations (NSE) for applications in nonhydrostatic atmospheric modeling. The compressible NSE in nonhydrostatic atmospheric modeling include buoyancy terms that require special handling if one wishes to extract the Schur complement form of the linear implicit problem. We present results for five different forms of the compressible NSE and describe in detail how to formulate the semi-implicit time-integration method for these equations. Finally, we compare all five equations and compare the semi-implicit formulations of these equations both using the Schur and No Schur forms against an explicit Runge-Kutta method. Our simulations show that, if efficiency is the main criterion, it matters which form of the governing equations you choose. Furthermore, the semi-implicit formulations are faster than the explicit Runge-Kutta method for all the tests studied, especially if the Schur form is used. While we have used the spectral element method for discretizing the spatial operators, the semi-implicit formulations that we derive are directly applicable to all other numerical methods. We show results for our five semi-implicit models for a variety of problems of interest in nonhydrostatic atmospheric modeling, including inertia-gravity waves, density current (i.e., Kelvin-Helmholtz instabilities), and mountain test cases; the latter test case requires the implementation of nonreflecting boundary conditions. Therefore, we show results for all five semi-implicit models using the appropriate boundary conditions required in nonhydrostatic atmospheric modeling: no-flux (reflecting) and nonreflecting boundary conditions (NRBCs). It is shown that the NRBCs exert a strong impact on the accuracy and efficiency of the models.

Analysis of adaptive mesh refinement for IMEX discontinuous Galerkin solutions of the compressible Euler equations with application to atmospheric simulations

The resolutions of interests in atmospheric simulations require prohibitively large computational resources. Adaptive mesh refinement (AMR) tries to mitigate this problem by putting high resolution in crucial areas of the domain. We investigate the performance of a tree-based AMR algorithm for the high order discontinuous Galerkin method on quadrilateral grids with non-conforming elements. We perform a detailed analysis of the cost of AMR by comparing this to uniform reference simulations of two standard atmospheric test cases: density current and rising thermal bubble. The analysis shows up to 15 times speed-up of the AMR simulations with the cost of mesh adaptation below 1% of the total runtime. We pay particular attention to the implicit-explicit (IMEX) time integration methods and show that the ARK2 method is more robust with respect to dynamically adapting meshes than BDF2. Preliminary analysis of preconditioning reveals that it can be an important factor in the AMR overhead. The compiler optimizations provide significant runtime reduction and positively affect the effectiveness of AMR allowing for speed-ups greater than it would follow from the simple performance model.