A discontinuous/continuous low order finite element shallow water model on the sphere
Introduction
Finite element (FE) schemes are successfully employed in numerous computational fluid dynamics applications. FE methods rely on a powerful mathematical apparatus and offer a good representation of the physical fields that are approximated using sets of basis functions of chosen accuracy. FE methods are applicable to unstructured grids and can accommodate grid refinement.
Today most dynamical cores of global atmosphere and ocean models are based on finite difference, or spectral transform schemes – mostly due to the ease of implementation and computational efficiency [4], [26]. Compared to those discretization methods, FE schemes allow more flexibility with regard to unstructured or locally refined grids and can approximate physical fields with higher order polynomials. There are dynamical cores of global ocean or atmosphere models based on FE schemes using continuous low order finite elements methods [10], spectral elements [33], or discontinuous Galerkin (DG) FE methods [11], [13], [23], [19].
Global scale presents special challenges for the dynamical core of an atmospheric or ocean model. In a finite element framework, four important properties need to be considered:
- 1.
The core needs to enable high resolution simulations.
- 2.
To support a smooth representation of topography and coastlines, one should be able to use very small grid cells. The second property is different to the first one since the internal resolution inside the grid cells can vary for different discretization schemes.
- 3.
The core must be able to represent the geostrophic balance between the Coriolis force and the pressure gradient force.
- 4.
The element must satisfy the Ladyzhenskaya–Babuska–Brezzi-condition (LBB-condition) – a necessary condition for stability of the discretization scheme.
The first two properties suggest low order elements. Unfortunately, it is difficult to find a low order element that fulfills the properties three and four at the same time. Failing in one of them leads to spurious modes [21]. Therefore, continuous low order finite element methods often require stabilization schemes [10].
A hybrid FE approach that combines a second order continuous representation for the scalar field with a first order discontinuous representation for the velocity field was proposed as a potential candidate to combine the four properties above. It was shown for the linear equations that the element is LBB-stable and is able to represent the geostrophic balance [6], [7]. The analysis of the element in terms of a Helmholtz decomposition for the linearized shallow-water equations showed that the element has no spurious pressure or Rossby modes. Spurious modes occur in the least harmful place via inertia oscillations that do not propagate [5].
In this paper we contribute to the understanding of the element. We use the element to solve the spherical shallow water equations with and without viscous dissipation and on unbounded domains as well as on domains with lateral boundaries. The solutions are analyzed with respect to their convergence, conservation properties, and energy spectra.
While there exists extensive literature on how to obtain stable finite element models when only continuous or only discontinuous field representations are used, it is not obvious how to obtain an adequate model setup for the hybrid continuous/discontinuous -element. As opposed to continuous FE configurations, we need to solve a Riemann problem at cell boundaries due to the discontinuous velocity field, and, for a stable model setup, we need to perform a partial integration of the nonlinear divergence term in the scalar equation. As opposed to typical DG approaches, we use the non-conservative form of the shallow water equations and non-orthogonal Lagrange polynomials as basis functions. Since the continuous scalar field is well defined on cell boundaries, our Riemann problem reduces to two dimensions for the two components of the velocity field, this is different to DG methods. Except for the Riemann solver, our FE configuration is similar to the one presented in [9] for a FE model on the plane. In [9] the weak form of the equations was derived from the non-conservative shallow water equations while the three-dimensional Riemann solver for velocity and scalar fields was deduced from the conservative form.
The introduction of the spherical geometry for global FE methods requires some care. To develop the physical fields into sets of basis functions, the triangles in physical space need to be mapped onto a reference triangle on which the basis functions are defined. On the sphere, the geometry of the physical triangles is given by trigonometric functions and cannot be mapped exactly onto a reference triangle, which is typically defined on the plane. In the literature, three different approaches can be found to introduce curved manifolds, such as the sphere, to finite element models. In the first approach, the differential equations are written for the curved manifold. This is done in the case of the cubed sphere, for example [28]. In the second approach, the differential equations are formulated in the three-dimensional space and a global Cartesian coordinate system is considered. The vector fields are forced to stay on the manifold via constraints [12]. In the third approach, the vector fields are written in the local tangent basis while the fluxes and spatial operators are expressed in the three dimensional Cartesian basis [3].
We propose the use of the stereographic projection to introduce the spherical geometry to global finite element models. The projection is part of the first approach mentioned above. It has already been used for finite difference and finite volume methods ([25], [8], [17], CCSR Ocean Component Model). In the stereographic projection, the sphere is projected from one of the poles onto a plane at the opposite side of the sphere. The pole itself is mapped to infinity. Global circulation models either use a combination of two stereographic projections from each pole that are connected at the equator or two stereographic projection from each pole that are connected to a Mercator projection in tropical regions [2], [25].
For an ocean model, it is sufficient to use one projection since the pole can be placed on land. For a global atmosphere model, we use two projections from each pole that are connected at the equator. The scaling factors that appear can be represented exactly by second order polynomials. The projection leads to a manageable number of additional terms in the differential equations.
In section two, we give an overview of the model setup, including the equations of motion, the discontinuous/continuous finite element discretization, and the incorporation of the spherical geometry. In section three, we apply the model to the standard test set for global shallow water models and a Stommel gyre ocean test case. In section four, we present a summary and conclusions.
Section snippets
The viscous shallow water equations
The viscous shallow water equations can be written in several different forms which are analytically equivalent but have different properties when discretized.
One of those forms is the non-conservative formwhere u is the two dimensional velocity vector, f is the Coriolis parameter, k is the vertical unit vector, g is the gravitational acceleration, ν is the eddy viscosity, τ is a forcing term (for example, bottom friction or wind stress in ocean
Numerical results
We evaluate our model using the standard suite of global shallow water test cases for the inviscid equation (ν = 0) as well as a non-linear Stommel gyre test with boundaries for the viscous case.
Geodesic grids
For our numerical experiments, we use icosahedral geodesic grids [1], due to their quasiuniform coverage of the sphere. In principle, each spherical triangular grid could be used for simulations. The grids can be refined to the next level by bisecting the sides of the triangles; the new
Conclusions
We investigate the potential applicability of the finite element pair in atmosphere and ocean modeling. We present a stable discretization for the viscous and inviscid shallow water equations on the sphere. We do not observe any spurious modes that necessitate additional diffusion or stabilization schemes.
Regarding the convergence properties and kinetic energy spectra, our discretization shows the expected behavior. A comparison to the ICON shallow water model shows very promising
Acknowledgments
We thank Hui Wan and the two anonymous reviewers for their useful and constructive comments, which helped to improve earlier versions of the paper.
References (35)
- et al.
High-order discontinuous Galerkin schemes on general 2d manifolds applied to the shallow water equations
Journal of Computational Physics
(2009) A numerical method for the study of the circulation of the world ocean
Journal of Computational Physics
(1969)- et al.
Numerical wave propagation for the triangular p1dg-p2 finite element pair
Journal of Computational Physics
(2011) - et al.
A mixed discontinuous/continuous finite element pair for shallow-water ocean modelling
Ocean Modelling
(2009) - et al.
Lbb stability of a mixed Galerkin finite element pair for fluid flow simulations
Journal of Computational Physics
(2009) - et al.
A finite-element ocean model: principles and evaluation
Ocean Modelling
(2004) - et al.
Nodal high-order discontinuous Galerkin methods for the spherical shallow water equations
J. Comput. Phys.
(2002) Lagrange-Galerkin methods on spherical geodesic grids
Journal of Computational Physics
(1997)High-order triangle-based discontinuous Galerkin methods for hyperbolic equations on a rotating sphere
Journal of Computational Physics
(2006)- et al.
Spatial discretization of the shallow water equations in spherical geometry using Osher’s scheme
Journal of Computational Physics
(2000)
A discontinuous Galerkin method for the shallow water equations in spherical triangular coordinates
Journal of Computational Physics
Unsteady analytical solutions of the spherical shallow water equations
Journal of Computational Physics
The ‘cubed sphere’: a new method for the solution of partial differential equations in spherical geometry
Journal of Computational Physics
An explicit expression for the penalty parameter of the interior penalty method
Journal of Computational Physics
The spectral element method for the shallow water equations on the sphere
Journal of Computational Physics
High-order finite-volume methods for the shallow-water equations on the sphere
Journal of Computational and Physics
A standard test set for numerical approximations to the shallow water equations in spherical geometry
Journal of Computational and Physics
Cited by (12)
A conservative discretization of the shallow-water equations on triangular grids
2018, Journal of Computational PhysicsNumerical Methods for the Nonlinear Shallow Water Equations
2017, Handbook of Numerical AnalysisCitation Excerpt :To remove its disadvantage of pole-singularity, other choices of meshes have been proposed and studied, including the Yin–Yang mesh, cubed-sphere mesh, among others. Well-balanced and positivity-preserving methods under these meshes, including spectral element, finite volume, continuous low order finite elements and DG methods, have also been investigated, and we refer to Taylor et al. (1997), Giraldo et al. (2002), Nair et al. (2005), Verkley (2009), Ullrich et al. (2010) and Duben et al. (2012) for recent development on numerical methods for the SWEs on the sphere. Solving the two-layer SWEs is a challenging problem due to several reasons: they are only conditionally hyperbolic; they contain nonconservative product terms; they admit steady state to be exactly preserved; their water heights should stay nonnegative, and their eigenstructure cannot be obtained in explicit form.
1-D coupled non-equilibrium sediment transport modeling for unsteady flows in the discontinuous Galerkin framework
2016, Journal of HydrodynamicsA novel Local Time Stepping algorithm for shallow water flow simulation in the discontinuous Galerkin framework
2016, Applied Mathematical ModellingCitation Excerpt :Shallow water flow equations, also known as Saint–Venant equations, are widely used in hydraulic and river engineering problems. Many researchers, with the use of different numerical methods such as finite volume, finite difference, finite element, hybrid finite element/volume, smoothed-particle hydrodynamics, and lattice Boltzmann, have successfully solved these equations [9–15]. In recent years, the discontinuous Galerkin (DG) finite element method has gained popularity in numerical modeling of shallow water flows [16–26].