A stable discontinuous Galerkin method for linear elastodynamics in 3D geometrically complex media using physics based numerical fluxes

Abstract

High order accurate and explicit time-stable solvers are well suited for hyperbolic wave propagation problems.For the complexities of real geometries, internal interfaces, nonlinear boundary and interface conditions, discontinuities and sharp wave fronts become fundamental features of the solutions. These are also effects of the presence of disparate spatial and temporal scales, present in real media and sources. As a result high order accuracy, geometrically flexible and adaptive numerical algorithms are critical for high fidelity and efficient simulations of wave phenomena in many applications.Using a physics-based numerical penalty-flux, we develop a provably energy-stable discontinuous Galerkin approximation of the elastic wave equation in complex and discontinuous media. By construction, our numerical flux is upwind and yields a discrete energy estimate analogous to the continuous energy estimate.The discrete energy estimate holds for conforming and non-conforming curvilinear elements. The ability to handle non-conforming curvilinear meshes allows for flexible adaptive mesh refinement strategies. The numerical scheme has been implemented in ExaHyPE, a simulation engine for hyperbolic PDEs on adaptive Cartesian meshes, for exascale supercomputers.We present 3D numerical experiments demonstrating stability and high order accuracy. Finally, we present a regional geophysical wave propagation problem in a 3D Earth model with geometrically complex free-surface topography.

Vectorization and Minimization of Memory Footprint for Linear High-Order Discontinuous Galerkin Schemes

“…We evaluate the performance of all the above-described kernel optimizations when simulating the elastic wave equations in first order formulation on curvilinear boundary-fitted meshes, as described in [8]. The equations are characterized by three quantities for particle velocity and six variables for the stress tensor.…”

“…ExaHyPE's performance is heavily dominated by the Space Time Predictor (STP), which computes an entirely elementlocal time extrapolation of the solution. For linear PDE systems, such as in the context of seismic simulations [8], the STP is computed via a Cauchy-Kowalewsky scheme, which requires tensor operations that imply calls to PDE-specific user functions. This generates conflicting requirements on the ExaHyPE API: the data layout needed for optimization of the tensor operations and that needed for vectorization of the user functions differs (AoS vs. SoA).…”

A stable discontinuous Galerkin method for the perfectly matched layer for elastodynamics in first order form

“…Computational strategies based on the discontinuous Galerkin method (DG method) are desirable for large scale numerical simulation of wave phenomena occurring in many applications [17,33,35,38]. However, one of the main features of propagating waves is that they can propagate long distances relative to their characteristic dimension, the wavelength.…”

“…To begin, we develop a DG numerical method for the linear elastodynamics equation using physically motivated numerical flux and penalty parameters, which are compatible with all well-posed, internal and external, boundary conditions [24,17]. When the PML damping vanishes, by construction, our choice of penalty parameters yield an upwind scheme and a discrete energy estimate analogous to the continuous energy estimate.…”

“…We discretize in time using the high order arbitrary derivative (ADER) time stepping scheme [17,32,31] of the same order of accuracy with the spatial discretization. By combining the DG spatial approximation with the high order ADER time stepping scheme and the accuracy of the PML we obtain an arbitrarily accurate wave propagation solver in the time domain.…”