High Order Discontinuous Galerkin Schemes Research Articles

Implementation of the transition model for high order discontinuous Galerkin method with hybrid discretization strategy

The γ-Reθt⌢ transition model has been applied within the framework of high order discontinuous Galerkin (DG) discretization of Reynolds-averaged Navier-Stokes (RANS) equations. A hybrid discretization strategy that takes advantage of all available information from DG method is suggested for the implementation of the transition model. Some techniques for the spatial and temporal discretization are utilized to enhance the effectiveness of the method. The resulting approaches achieve high accuracy of transition prediction and meanwhile do not increase notable computational challenge when compared to DG discretization of fully coupled RANS/transition system. We carry out both steady and unsteady transitional simulations to examine the performance of the designed method. Numerical results have demonstrated that the present approaches are capable of capturing the transition flow characteristics as well as improving transition prediction accuracy with high order DG schemes.

An efficient sliding mesh interface method for high-order discontinuous Galerkin schemes

Sliding meshes are a powerful method to treat deformed domains in computational fluid dynamics, where different parts of the domain are in relative motion. In this paper, we present an efficient implementation of a sliding mesh method into a discontinuous Galerkin compressible Navier-Stokes solver and its application to a large eddy simulation of a 1-1/2 stage turbine. The method is based on the mortar method and is high-order accurate. It can handle three-dimensional sliding mesh interfaces with various interface shapes. For plane interfaces, which are the most common case, conservativity and free-stream preservation are ensured. We put an emphasis on efficient parallel implementation. Our implementation generates little computational and storage overhead. Inter-node communication via MPI in a dynamically changing mesh topology is reduced to a bare minimum by ensuring a priori information about communication partners and data sorting. We provide performance and scaling results showing the capability of the implementation strategy. Apart from analytical validation computations and convergence results, we present a wall-resolved implicit LES of the 1-1/2 stage Aachen turbine test case as a large scale practical application example.

Numerical study of two models for viscous compressible fluid flows

Our aim is to evaluate the modelling capability of a recently proposed (Eulerian) model for viscous and heat conducting compressible fluids. The new model differs from the Navier-Stokes in the structure of the diffusive terms. Most notably, it includes a mass diffusive mechanism. Despite the lack of well-posedness theory, one has to recognize that the standard Navier-Stokes equations are known to give reasonably accurate results for numerous standard problems. Hence, the new model can not deviate too much from the Navier-Stokes equations for such cases, is it to be an accurate model. The purpose of the paper is to compare the new model with Navier-Stokes and determine an unknown model parameter (α) in the new model such that the best matches are obtained. Using a high-order discontinuous Galerkin scheme, with careful implementations of the boundary conditions, we perform accurate simulations of the two models for a suite of test cases. For the grid converged results, we have, for each problem, compared the two models with respect to a number of measures such as cL, cD, cM, skin friction coefficients, pressure coefficients and Mach number or pressure isolines. The numerical results reveal a remarkable sameness (the differences are typically 1%) of both models when α=1.

A weight-adjusted discontinuous Galerkin method for wave propagation in coupled elastic-acoustic media

This paper presents a high-order discontinuous Galerkin (DG) scheme for the simulation of wave propagation through coupled elastic-acoustic media. We use a first-order stress-velocity formulation, and derive a simple upwind-like numerical flux which weakly imposes continuity of the normal velocity and traction at elastic-acoustic interfaces. When combined with easily invertible weight-adjusted mass matrices [1–3], the resulting method is efficient, consistent, and energy stable on curvilinear meshes and for arbitrary heterogeneous media, including anisotropy and sub-cell (micro) heterogeneities. We numerically verify the high order accuracy and stability of the proposed method, and investigate its performance for applications in photoacoustic tomography.

A high-order discontinuous Galerkin method for nonlinear sound waves

We propose a high-order discontinuous Galerkin scheme for nonlinear acoustic waves on polytopic meshes. To model sound propagation with losses through homogeneous media, we use Westervelt's nonlinear wave equation with strong damping. Challenges in the numerical analysis lie in handling the nonlinearity in the model, which involves the derivatives in time of the acoustic velocity potential, and in preventing the equation from degenerating. We rely in our approach on the Banach fixed-point theorem combined with a stability and convergence analysis of a linear wave equation with a variable coefficient in front of the second time derivative. By doing so, we derive an a priori error estimate for Westervelt's equation in a suitable energy norm for the polynomial degree p≥2. Numerical experiments carried out in two-dimensional settings illustrate the theoretical convergence results. In addition, we demonstrate efficiency of the method in a three-dimensional domain with varying medium parameters, where we use the discontinuous Galerkin approach in a hybrid way.

Continuum electromagnetic gyrokinetic simulations of turbulence in the tokamak scrape-off layer and laboratory devices

We present algorithms and results from Gkeyll, a full-f continuum, electromagnetic gyrokinetic code, designed to study turbulence in the edge region of fusion devices. The edge is computationally very challenging, requiring robust algorithms that can handle large-amplitude fluctuations and stable interactions with plasma sheaths. We present an energy-conserving high-order discontinuous Galerkin scheme that solves gyrokinetic equations in Hamiltonian form. Efficiency is improved by a careful choice of basis functions and automatically generated computation kernels. Previous verification tests were performed in the straight-field-line large plasma device [Shi et al., J. Plasma Phys. 83, 905830304 (2017)] and the Texas Helimak, a simple magnetized torus [Bernard et al., Phys. Plasmas 26, 042301 (2019)], including the effect of end-plate biasing on turbulence. Results for the scrape-off layer for NSTX parameters with a model helical magnetic geometry with bad curvature have been obtained [Shi et al., Phys. Plasmas 26, 012307 (2019)]. In this paper, we present algorithms for the two formulations of electromagnetic gyrokinetics: the Hamiltonian and the symplectic. We describe each formulation and show results of benchmark tests. Although our scheme works for the Hamiltonian formulation, the presence of spurious numerical modes for high-β and large k⊥2ρs2 regimes shows that the symplectic formulation is more robust. We then review our recent algorithm for the symplectic formulation [Mandell et al., J. Plasma Phys. 86, 905860109 (2020)], along with example application of this new capability. Maintaining positivity of the distribution function can be challenging, and we describe a new and novel exponential recovery based algorithm to address this.

Controlling oscillations in high-order Discontinuous Galerkin schemes using artificial viscosity tuned by neural networks

High-order numerical solvers for conservation laws suffer from Gibbs phenomenon close to discontinuities, leading to spurious oscillations and a detrimental effect on the solution accuracy. A possible strategy to reduce it comprises adding a suitable amount of artificial dissipation. Although several viscosity models have been proposed in the literature, their dependence on problem-dependent parameters often limits their performances. Motivated by the objective to construct a universal artificial viscosity method, we propose a new technique based on neural networks, integrated into a Runge-Kutta Discontinuous Galerkin solver. Numerical results are presented to demonstrate the performance of this network-based technique. We show that it is able both to guarantee optimal accuracy for smooth problems, and to accurately detect discontinuities, where dissipation has to be injected. A comparison with some classical models is carried out, showing that the network-based model is always at par with the best among the traditional optimized models, independently of the selected problem and parameters.

Entropy stable numerical approximations for the isothermal and polytropic Euler equations

In this work we analyze the entropic properties of the Euler equations when the system is closed with the assumption of a polytropic gas. In this case, the pressure solely depends upon the density of the fluid and the energy equation is not necessary anymore as the mass conservation and momentum conservation then form a closed system. Further, the total energy acts as a convex mathematical entropy function for the polytropic Euler equations. The polytropic equation of state gives the pressure as a scaled power law of the density in terms of the adiabatic index gamma . As such, there are important limiting cases contained within the polytropic model like the isothermal Euler equations (gamma {=}1) and the shallow water equations (gamma {=}2). We first mimic the continuous entropy analysis on the discrete level in a finite volume context to get special numerical flux functions. Next, these numerical fluxes are incorporated into a particular discontinuous Galerkin (DG) spectral element framework where derivatives are approximated with summation-by-parts operators. This guarantees a high-order accurate DG numerical approximation to the polytropic Euler equations that is also consistent to its auxiliary total energy behavior. Numerical examples are provided to verify the theoretical derivations, i.e., the entropic properties of the high order DG scheme.

A high-order conservative remap for discontinuous Galerkin schemes on curvilinear polygonal meshes

A data transfer (called later remap) of physical fields between two meshes is an important step of arbitrary Lagrangian-Eulerian (ALE) simulations. This step is challenging for high-order discontinuous Galerkin schemes since the Lagrangian flow motion leads to high-order meshes with curved faces. It becomes even more challenging for unstructured polygonal meshes that do not have a polynomial map from the reference to a current cell. We propose and analyze a new framework to create remap schemes on curvilinear polygonal meshes based on the theory of virtual element projectors. We derive a conservative remap scheme that is high-order accurate in space and time. The properties of this scheme are studied numerically for smooth and discontinuous fields on unstructured quadrilateral and polygonal meshes.

Utilization of the Brinkman Penalization to Represent Geometries in a High-Order Discontinuous Galerkin Scheme on Octree Meshes

We investigate the suitability of the Brinkman penalization method in the context of a high-order discontinuous Galerkin scheme to represent wall boundaries in compressible flow simulations. To evaluate the accuracy of the wall model in the numerical scheme, we use setups with symmetric reflections at the wall. High-order approximations are attractive as they require few degrees of freedom to represent smooth solutions. Low memory requirements are an essential property on modern computing systems with limited memory bandwidth and capability. The high-order discretization is especially useful to represent long traveling waves, due to their small dissipation and dispersion errors. An application where this is important is the direct simulation of aeroacoustic phenomena arising from the fluid motion around obstacles. A significant problem for high-order methods is the proper definition of wall boundary conditions. The description of surfaces needs to match the discretization scheme. One option to achieve a high-order boundary description is to deform elements at the boundary into curved elements. However, creating such curved elements is delicate and prone to numerical instabilities. Immersed boundaries offer an alternative that does not require a modification of the mesh. The Brinkman penalization is such a scheme that allows us to maintain cubical elements and thereby the utilization of efficient numerical algorithms exploiting symmetry properties of the multi-dimensional basis functions. We explain the Brinkman penalization method and its application in our open-source implementation of the discontinuous Galerkin scheme, Ateles. The core of this presentation is the investigation of various penalization parameters. While we investigate the fundamental properties with one-dimensional setups, a two-dimensional reflection of an acoustic pulse at a cylinder shows how the presented method can accurately represent curved walls and maintains the symmetry of the resulting wave patterns.

An efficient hyperbolic relaxation system for dispersive non-hydrostatic water waves and its solution with high order discontinuous Galerkin schemes

In this paper we propose a novel set of first-order hyperbolic equations that can model dispersive non-hydrostatic free surface flows. The governing PDE system is obtained via a hyperbolic approximation of the family of non-hydrostatic free-surface flow models recently derived by Sainte-Marie et al. in [1]. Our new hyperbolic reformulation is based on an augmented system in which the divergence constraint of the velocity is coupled with the other conservation laws via an evolution equation for the depth-averaged non-hydrostatic pressure, similar to the hyperbolic divergence cleaning applied in generalized Lagrangian multiplier methods (GLM) for magnetohydrodynamics (MHD). We suggest a formulation in which the divergence errors of the velocity field are transported with a large but finite wave speed that is directly related to the maximal eigenvalue of the governing PDE.We then use arbitrary high order accurate (ADER) discontinuous Galerkin (DG) finite element schemes with an a posteriori subcell finite volume limiter in order to solve the proposed PDE system numerically. The final scheme is highly accurate in smooth regions of the flow and very robust and positive preserving for emerging topographies and wet-dry fronts. It is well-balanced making use of a path-conservative formulation of HLL-type Riemann solvers based on the straight line segment path. Furthermore, the proposed ADER-DG scheme with a posteriori subcell finite volume limiter adapts very well to modern GPU architectures, resulting in a very accurate, robust and computationally efficient computational method for non-hydrostatic free surface flows. The new model proposed in this paper has been applied to idealized academic benchmarks such as the propagation of solitary waves, as well as to more challenging physical situations that involve wave runup on a shore including wave breaking in both one and two space dimensions. In all cases the achieved agreement with analytical solutions or experimental data is very good, thus showing the validity of both, the proposed mathematical model and the numerical solution algorithm.

A posteriori correction of high-order discontinuous Galerkin scheme through subcell finite volume formulation and flux reconstruction

In this paper, we present a new limiter for discontinuous Galerkin (DG) schemes, based on subcell resolution through reconstructed flux correction, for hyperbolic conservation laws. This limitation technique is constructed by means of a subcell Finite Volume (FV) formulation. The algorithm is thus simple, and is able to retain the very high accurate subcell resolution of DG schemes.The main idea of this correction procedure is to preserve as much as possible the high accuracy and the very precise subcell resolution of DG schemes. Consequently, an a posteriori correction will only be applied locally at the subcell scale where it is needed, but still ensuring the scheme conservation. To do so, we first reformulate DG schemes as a subcell FV scheme provided the use of the correct numerical flux. This very simple development provides us with the so-called DG reconstructed flux. This theoretical part will serve as a basis for our limiter framework.Practically, at each time step, we compute a DG candidate solution and check if this solution is admissible (for instance positive, non-oscillating, …). If it is the case, we go further in time. Otherwise, we return to the previous time step and correct locally, at the subcell scale, the numerical solution. This is why it is refereed to as a posteriori limitation. To this end, each cell is subdivided into subcells. Then, if the solution is locally detected as bad, we substitute the DG reconstructed flux on the subcell boundaries by a robust first-order or second-order TVD numerical flux. And for subcell detected as admissible, we keep the high-order reconstructed flux which allows us to retain the very high accurate resolution and conservation of the DG scheme. Furthermore, only the solution inside troubled subcells and its first neighbors will have to be recomputed, elsewhere the solution remains unchanged. Numerical results on various type problems and test cases will be presented, both in 1D and 2D on Cartesian grids, to assess the very good performance of the design limiting algorithm.

Mesh Optimization for Discontinuous Galerkin Methods Using a Continuous Mesh Model

A method for anisotropic mesh adaptation and optimization for high-order discontinuous Galerkin schemes is presented. Given the total number of degrees of freedom, a metric-based method is proposed, which aims to globally optimize the mesh with respect to the norm of the error. This is done by minimizing a suitable error model associated with the approximation space. Advantages of using a metric-based method in this context are several. First, it facilitates changing and manipulating the mesh in a general anisotropic way. Second, defining a suitable continuous interpolation operator allows the use of an analytic optimization framework that operates on the metric field, rather than the discrete mesh. The formulation of the method is presented as well as numerical experiments in the context of convection–diffusion systems.

An invariant‐region‐preserving limiter for DG schemes to isentropic Euler equations

In this article, we introduce an invariant‐region‐preserving (IRP) limiter for the p‐system and the corresponding viscous p‐system, both of which share the same invariant region. Rigorous analysis is presented to show that for smooth solutions the order of approximation accuracy is not destroyed by the IRP limiter, provided the cell average stays away from the boundary of invariant region. Moreover, this limiter is explicit, and easy for computer implementation. A generic algorithm incorporating the IRP limiter is presented for high order finite volume type schemes as long as the evolved cell average of the underlying scheme stays strictly within the invariant region. For high order discontinuous Galerkin (DG) schemes to the p‐system, sufficient conditions are obtained for cell averages to stay in the invariant region. For the viscous p‐system, we design both second and third order IRP DG schemes. Numerical experiments are provided to test the proven properties of the IRP limiter and the performance of IRP DG schemes.

Conformal and covariant Z4 formulation of the Einstein equations: Strongly hyperbolic first-order reduction and solution with discontinuous Galerkin schemes

We present a strongly hyperbolic first-order formulation of the Einstein equations based on the conformal and covariant Z4 system (CCZ4) with constraint-violation damping, which we refer to as FO-CCZ4. As CCZ4, this formulation combines the advantages of a conformal and traceless formulation, with the suppression of constraint violations given by the damping terms, but being first order in time and space, it is particularly suited for a discontinuous Galerkin (DG) implementation. The strongly hyperbolic first-order formulation has been obtained by making careful use of first and second-order ordering constraints. A proof of strong hyperbolicity is given for a selected choice of standard gauges via an analytical computation of the entire eigenstructure of the FO-CCZ4 system. The resulting governing partial differential equations system is written in non-conservative form and requires the evolution of 58 unknowns. A key feature of our formulation is that the first-order CCZ4 system decouples into a set of pure ordinary differential equations and a reduced hyperbolic system of partial differential equations that contains only linearly degenerate fields. We implement FO-CCZ4 in a high-order path-conservative arbitrary-high-order-method-using-derivatives (ADER)-DG scheme with adaptive mesh refinement and local time-stepping, supplemented with a third-order ADER-WENO subcell finite-volume limiter in order to deal with singularities arising with black holes. We validate the correctness of the formulation through a series of standard tests in vacuum, performed in one, two and three spatial dimensions, and also present preliminary results on the evolution of binary black-hole systems. To the best of our knowledge, these are the first successful three-dimensional simulations of moving punctures carried out with high-order DG schemes using a first-order formulation of the Einstein equations.

A positivity-preserving high order discontinuous Galerkin scheme for convection–diffusion equations

For constructing high order accurate positivity-preserving schemes for convection–diffusion equations, we construct a simple positivity-preserving diffusion flux. Discontinuous Galerkin (DG) schemes with such a positivity-preserving diffusion flux are nonlinear schemes, which can be regarded as a reduction of the high order positivity-preserving DG schemes for compressible Navier–Stokes equations in [1] to scalar diffusion operators. In this paper we focus on the local DG method to discuss how to apply such a flux. A limiter on the auxiliary variable for approximating the gradient of the solution must be used so that the diffusion flux is positivity-preserving in the sense that DG schemes with this flux satisfies a weak positivity property. Together with a positivity-preserving limiter, high order DG schemes with strong stability preserving time discretizations can be rendered positivity-preserving without losing conservation or high order accuracy for convection–diffusion problems with periodic boundary conditions or a special class of Dirichlet or Neumann boundary conditions. Numerical tests on a few parabolic equations and an application to modeling electrical discharges are shown to demonstrate the performance of this scheme.

ADER discontinuous Galerkin schemes for general-relativistic ideal magnetohydrodynamics

We present a new class of high-order accurate numerical algorithms for solving the equations of general-relativistic ideal magnetohydrodynamics in curved spacetimes. In this paper we assume the background spacetime to be given and static, i.e., we make use of the Cowling approximation. The governing partial differential equations are solved via a new family of fully-discrete and arbitrary high-order accurate path-conservative discontinuous Galerkin (DG) finite-element methods combined with adaptive mesh refinement and time accurate local timestepping. In order to deal with shock waves and other discontinuities, the highorder DG schemes are supplemented with a novel a-posteriori subcell finite-volume limiter, which makes the new algorithms as robust as classical second-order total-variation diminishing finite-volume methods at shocks and discontinuities, but also as accurate as unlimited high-order DG schemes in smooth regions of the flow. We show the advantages of this new approach by means of various classical two- and three-dimensional benchmark problems on fixed spacetimes. Finally, we present a performance and accuracy comparisons between Runge-Kutta DG schemes and ADER high-order finite-volume schemes, showing the higher efficiency of DG schemes.

Invariant-region-preserving DG methods for multi-dimensional hyperbolic conservation law systems, with an application to compressible Euler equations

An invariant-region-preserving (IRP) limiter for multi-dimensional hyperbolic conservation law systems is introduced, as long as the system admits a global invariant region which is a convex set in the phase space. It is shown that the order of approximation accuracy is not destroyed by the IRP limiter, provided the cell average is away from the boundary of the convex set. Moreover, this limiter is explicit, and easy for computer implementation. A generic algorithm incorporating the IRP limiter is presented for high order finite volume type schemes. For arbitrarily high order discontinuous Galerkin (DG) schemes to hyperbolic conservation law systems, sufficient conditions are obtained for cell averages to remain in the invariant region provided the projected one-dimensional system shares the same invariant region as the full multi-dimensional hyperbolic system does. The general results are then applied to both one and two dimensional compressible Euler equations so to obtain high order IRP DG schemes. Numerical experiments are provided to validate the proven properties of the IRP limiter and the performance of IRP DG schemes for compressible Euler equations.

Semi-implicit staggered discontinuous Galerkin schemes for axially symmetric viscous compressible flows in elastic tubes

We propose a novel family of staggered semi-implicit discontinuous Galerkin (DG) finite element schemes for the simulation of axially symmetric, weakly compressible and laminar viscous flows in elastic pipes. The equation of state (EOS) of the fluid is assumed to be barotropic and two different mathematical models derived from the compressible Navier-Stokes equations are considered in this paper.The first model describes cross-sectionally averaged 1D flows, including steady and frequency-dependent wall friction effects. The novelty of our numerical method compared to standard DG schemes consists in the use of a staggered mesh, where the pressure is defined over a primary grid and the velocity field is defined on edge-based staggered dual control volumes. This approach is well known from classical semi-implicit finite difference schemes for the incompressible Navier–Stokes equations, but it is still quite unusual for high order DG schemes. The continuity equation is integrated over the control volumes that belong to the main grid, while the momentum equation is integrated over the elements of the edge-based staggered dual grid. The nonlinear convective terms are discretized explicitly, while the pressure gradient and the mass flux are discretized implicitly. Up to second order of accuracy in time can be achieved with the so-called θ-method. Inserting the discrete momentum equation in the discrete mass conservation equation leads to a mildly nonlinear algebraic system for the degrees of freedom of the pressure. Such mildly nonlinear systems can be very efficiently solved using the Newton algorithm of Brugnano and Casulli. We observe that the linear part of the mildly nonlinear system is symmetric and positive definite.The second model is derived from the compressible Navier-Stokes equations in cylindrical coordinates. Assuming hydrostatic flow with constant pressure inside each cross section as well as axial symmetry, only the terms in the axial and the radial direction need to be considered. Therefore, we call the second model the 2Dxr model. Also in this case we use a staggered mesh for pressure and velocity and thus the same philosophy as for the 1D model can be applied to obtain the discrete pressure system. For the 2Dxr model a staggered DG scheme is also applied for the computation of the viscous stress tensor in the discrete momentum equation. However, in radial direction the resulting linear system for the friction terms is not symmetric and is thus solved using the Thomas algorithm for block three-diagonal systems.The use of a semi-implicit DG scheme leads to a very mild CFL condition based only on the fluid velocity and not on the sound speed, which makes the method very efficient, in particular in the limit cases when the speed of sound of the fluid tends to infinity (incompressible fluid) and in the rigid case where the wall strain of the pipe tends to zero. In addition, at every Newton step a symmetric positive definite and well conditioned block three-diagonal linear system is solved for the pressure, using a matrix-free conjugate gradient method. Moreover, when the polynomial degree of the basis and test functions is equal to zero the schemes reduce to classical semi-implicit finite volume methods.While in the 2Dxr model the viscous effects in radial direction are directly obtained from first principles via the Navier-Stokes equations, the 1D model requires an additional closure relation for the wall friction. For both models we perform several tests in order to validate the numerical methods for steady and unsteady flows of compressible and nearly incompressible fluids in elastic and rigid tubes. We also provide numerical convergence results in order to show that the developed schemes achieve high order of accuracy in space.

Efficient Explicit Time Stepping of High Order Discontinuous Galerkin Schemes for Waves

This work presents algorithms for the efficient implementation of discontinuous Galerkin methods with explicit time stepping for acoustic wave propagation on unstructured meshes of quadrilaterals or hexahedra. A crucial step towards efficiency is to evaluate operators in a matrix-free way with sum-factorization kernels. The method allows for general curved geometries and variable coefficients. Temporal discretization is carried out by low-storage explicit Runge-Kutta schemes and the arbitrary derivative (ADER) method. For ADER, we propose a flexible basis change approach that combines cheap face integrals with cell evaluation using collocated nodes and quadrature points. Additionally, a degree reduction for the optimized cell evaluation is presented to decrease the computational cost when evaluating higher order spatial derivatives as required in ADER time stepping. We analyze and compare the performance of state-of-the-art Runge-Kutta schemes and ADER time stepping with the proposed optimizations. ADER involves fewer operations and additionally reaches higher throughput by higher arithmetic intensities and hence decreases the required computational time significantly. Comparison of Runge-Kutta and ADER at their respective CFL stability limit renders ADER especially beneficial for higher orders when the Butcher barrier implies an overproportional amount of stages. Moreover, vector updates in explicit Runge--Kutta schemes are shown to take a substantial amount of the computational time due to their memory intensity.