A discontinuous Galerkin method for the shallow water equations in spherical triangular coordinates☆
Introduction
Modeling atmospheric flows for climate simulations as well as for weather prediction is a complex problem, due to the nonlinear structure of the dynamical and physical phenomena on widely varying spatial and temporal scales and their multi-scale interaction processes. Depending on the complexity of an atmospheric model the governing equations are the fundamental atmospheric conservation laws for mass, momentum, and energy or appropriate simplifications of them. If the regarded equation set is a hyperbolic system, energetic shocks can develop theoretically. Although this is usually not the case in atmospheric models, the discretization should represent regions of scale collapse and breaking waves generating discontinuities in the velocity field; the discrete conservation properties of the discontinuous Galerkin (DG) method are appropriate for this task.
The shallow water equations (SWE), valid for a homogeneous atmosphere with small vertical velocities and horizontal velocities independent in the vertical direction, constitute a hyperbolic system of conservation laws. It is one of the simplest nonlinear hyperbolic systems, covering important planetary atmospheric features, like the Rossby wave formation.
For the spherical SWE the spatial domain is the sphere S, a two-dimensional surface in . In a regional or mesoscale SWE model the momentum is a two-dimensional vector. In contrast, the Cartesian formulation of the spherical case in [43], [8], [44] represents the tangential momentum of the flow as a three-dimensional vector and includes a Lagrangian multiplier to constrain the momentum to be tangential. Applying this form to a numerical model usually leads to three momentum equations and requires a correction step to satisfy the constraint discretely, see e.g. [17]. Models in standard spherical coordinates satisfy a two-dimensional momentum representation but have to pay additional attention to phenomena near the poles due to singularities of the coordinate mapping, see e.g. [27].
The idea to avoid the drawbacks of the Cartesian and the spherical coordinates formulation is to represent the spherical SWE in terms of local coordinate transformations. On a cubed-sphere grid, a spherical quadrilateral grid, the following works [37], [33], [31], [45], [35] achieved a two-dimensional momentum representation avoiding any pole problem. Our new model achieves the same flexibility but on unstructured triangular grids using spherical triangular coordinates.
Numerous models on spherical triangles have been proposed in the last three decades. For example, the early work by Sadourny at al. [38] and Williamson [49] introduced to the atmospheric community the use of triangular grids based on the icosahedron to develop the underlying grid for the construction of finite difference operators. Work on triangular grids based on the icosahedron for discretizing the sphere lay dormant for another 15 years until the work was resumed by Baumgardner and Frederickson [1]. Ten years later work on these grids was resurrected by Heikes and Randall [20] and followed by Giraldo [13], Thuburn [47], Stuhne and Peltier [42], Tomita et al. [48], and Heinze and Hense [21]. All of these models rely directly on either the triangular grid being derived from the icosahedron or on a linear representation of the discrete operators; this way there is an easily computable dual grid which is based on hexagons (icosahedron) or the operators can be constructed using the vertices of an element which are co-planar (linear representation). However, for high-order operators on general triangular grids, one needs to construct the discrete spatial operators directly on the curved manifold which then requires the derivation of the Christoffel symbols from differential geometry (see [35] for a summary of the use of differential geometry for atmospheric flow).
One of the contributions of this manuscript is to show how to use these ideas for constructing high-order spatial operators on general triangulations on the sphere. While these ideas have been used extensively for quadrilateral-based grids, they have not been used at all (or not often); to our knowledge, there are currently no existing models on the sphere which use high-order discretizations on the triangle with the exception being those developed by the present authors. Läuter et al. [25] developed a shallow water model on the sphere using second order finite elements on dynamically adaptive triangular grids, while Giraldo and Warburton [17] and Giraldo [15] developed shallow water models on the sphere using up to 15 order spectral element and discontinuous Galerkin operators on unstructured triangular grids. Extracting the highlights of these three separate works results in a robust, accurate, and efficient model. However, in order to achieve this aim requires writing the equations directly on the parametric space defined by the curved spherical triangles. We consider the spherical SWE in flux form on the surface S using coordinate independent differential operators on S known from differential geometry. By using an appropriate local coordinate mapping on the curved spherical triangle E we can then obtain a two-dimensional representation of the tangential momentum vectors. This then allows us full grid independence such that our discrete operators can be constructed to arbitrarily high-order while doing so on generalized unstructured triangulations on the sphere (or any other curved manifold, for that matter). This approach allows us the same flexibility enjoyed with the Cartesian methods discussed previously while now only requiring tangential momentum equations.
Numerous numerical methods have been proposed for next generation global atmospheric models including finite volumes [27], [34], spectral elements [45], [11], [9], [17], and DG [16], [30], [15] methods. We have selected the DG method for our model because it allows us to achieve high-order accuracy as in spectral elements while conserving all quantities both locally and globally as in finite volumes, see the review in [5]. Furthermore, our use of unstructured triangular grids allows for much flexibility in future work on adaptivity.
The organization of this article is as follows. In Section 2 the governing spherical SWE are given using surface differential operators. Section 3 describes the numerical discretization by a Runge–Kutta discontinuous Galerkin method applying spherical triangular coordinates. In Section 4 the atmospheric model based on the discretization is validated in terms of standard tests from Williamson et al. [50], steady-state and unsteady analytical solutions and a barotropic instability generated by a small initial perturbation.
Section snippets
Spherical shallow water equations
The spherical SWE are a system of conservation laws for the geopotential layer depth (mass) and the flow momentum. Because the integration domain of the SWE is the sphere, a two-dimensional surface in , the system can be formulated in the surrounding Cartesian space , see Côté [8]. Côté’s formulation is equivalent to a conservative form of the SWE on the surface S, which is the formulation used further below.
Let us consider the sphere with the Earth’s radius , the
Discontinuous Galerkin method
The DG method is applied to the conservative form (1) of the SWE on the surface S. On each curved triangle (element) E of the grid tesselation spherical triangular coordinates are introduced, which are local coordinate mappings on E. The polynomial representation on each grid element uses high-order Lagrange polynomials based on specially chosen Lagrange points (see Section 3.2). This approach leads to the local representation of the tangential momentum fields by two components only. An
Numerical results
The RK-DG method described in Section 3 has been used to implement an atmospheric model which has been validated performing numerical experiments. The validation process is carried out in three steps. At first, standard model tests (flow over an isolated mountain and a Rossby–Haurwitz wave) from Williamson [50] are performed. After that, a convergence study considering steady-state and unsteady analytical solutions of the nonlinear SWE is done. Finally, a barotropic instability in a localized
Summary and outlook
A global model of the atmosphere has been developed based on the spherical shallow water equations and discretized by a Runge–Kutta discontinuous Galerkin method. The equations are formulated as a hyperbolic conservation law on the sphere, a two-dimensional surface in , using coordinate independent differential operators.
Spherical triangular coordinates, which are local coordinate mappings , and high-order polynomial spaces are defined on each curved element of a given spherical
Acknowledgements
The authors gratefully acknowledge the support of the Office of Naval Research Global through the Grant N00014-07-1-4038. The second author also gratefully acknowledges the support of the Office of Naval Research through PE-0602435N.
References (50)
- et al.
amatos: parallel adaptive mesh generator for atmospheric and oceanic simulation
Ocean Model.
(2005) - et al.
Interpolation theory over curved elements, with applications to finite element methods
Comput. Methods Appl. Mech. Engrg.
(1972) Monomial cubature rules since Stroud: a compilation – part 2
J. Comput. Appl. Math.
(1999)- et al.
Monomial cubature rules since Stroud: a compilation
J. Comput. Appl. Math.
(1993) - et al.
A semi-implicit discontinuous galerkin finite element method for the numerical solution of inviscid compressible flow
J. Comput. Phys.
(2004) Lagrange-Galerkin methods on spherical geodesic grids
J. Comput. Phys.
(1997)High-order triangle-based discontinuous Galerkin methods for hyperbolic equations on a rotating sphere
J. Comput. Phys.
(2006)- et al.
Nodal high-order discontinuous Galerkin methods for the spherical shallow water equations
J. Comput. Phys.
(2002) - et al.
A nodal triangle-based spectral element method for the shallow water equations on the sphere
J. Comput. Phys.
(2005) - et al.
Shock detection and limiting with discontinuous Galerkin methods for hyperbolic conservation laws
Appl. Numer. Math.
(2004)
Unsteady analytical solutions of the spherical shallow water equations
J. Comput. Phys.
A parallel adaptive barotropic model of the atmosphere
J. Comput. Phys.
The cubed sphere: a new method for the solution of partial differential equations in spherical geometry
J. Comput. Phys.
A wave propagation method for hyperbolic systems on the sphere
J. Comput. Phys.
A wave propagation algorithm for hyperbolic systems on curved manifolds
J. Comput. Phys.
Efficient implementation of essentially non-oscillatory shock-capturing schemes
J. Comput. Phys.
New icosahedral grid-point discretizations of the shallow water equations on the sphere
J. Comput. Phys.
The spectral element method for the shallow water equations on the sphere
J. Comput. Phys.
Shallow water model on a modified icosahedral geodesic grid by using spring dynamics
J. Comput. Phys.
A standard test set for numerical approximations to the shallow water equations in spherical geometry
J. Comput. Phys.
Icosahedral discretization of the two-sphere
SIAM J. Numer. Anal.
The Runge–Kutta local projection -discontinuous-Galerkin finite element method for scalar conservation laws
Math. Mod. Num. Anal.
Runge–Kutta discontinuous Galerkin methods for convection-dominated problems
J. Sci. Comput.
A Lagrange multiplier approach for the metric terms of semi-Lagrangian models on the sphere
Q.J.R. Meteorol. Soc.
High-resolution mesh convergence porperties and parallel efficiency of a spectral element atmospheric dynamical core
Int. J. High Perfor. Comput. Appl.
Cited by (47)
Reducing errors caused by geometrical inaccuracy to solve partial differential equations with moving frames on curvilinear domain
2022, Computer Methods in Applied Mechanics and EngineeringDiscontinuous Galerkin solver for the shallow-water equations in covariant form on the sphere and the ellipsoid
2020, Journal of Computational PhysicsIMEX HDG-DG: A coupled implicit hybridized discontinuous Galerkin and explicit discontinuous Galerkin approach for shallow water systems
2020, Journal of Computational PhysicsA conservative discretization of the shallow-water equations on triangular grids
2018, Journal of Computational PhysicsDispersion analysis of the P<inf>n</inf>−P<inf>n−1</inf><sup>DG</sup> mixed finite element pair for atmospheric modelling
2018, Journal of Computational PhysicsCitation Excerpt :However as noted by Adcroft et al. [1], Le Roux et al. [28] and others, the C-grid supports a stationary Coriolis mode due to the spatial averaging necessary to reconstruct the tangential velocity components which can become problematic for poorly resolved Rossby radius. In recent years high-order methods such as spectral elements [21,38], and discontinuous Galerkin (DG) based methods [18,20,32] have increased in popularity, due in part to the high level of computational efficiency obtainable on massively parallel supercomputers [32] as well as providing a high order of accuracy on non-uniform, non-orthogonal grids. This gives a greater flexibility in the choice of grids, thus allowing a more uniform tiling of the sphere or adaptive/variable resolution grids [24,7].
- ☆
This work relates to Department of the Navy Grant N00014-07-1-4038 issued by the Office of Naval Research Global. The United States Government has a royalty-free licence throughout the world in all copyrightable material contained herein.