Upwind summation-by-parts (SBP) finite difference methods for 3D seismic wave propagation in complex geometries

Seismic waves emanating from geophysical events propagate over hundreds to thousands of kilometers interacting with tectonic forces, geological structure, complicated topography and earthquake source processes on scales down to millimeters. Exploration seismology and natural earthquake hazard mitigation increasingly rely on multi-scale (0–20 Hz). High order accurate and computationally efficient numerical solvers become critical for efficient numerical simulation of strong ground motions from earthquakes.

Computational efficiency has continued to make the use of finite difference methods on structured grids attractive. An example is the petascale high order accurate SBP finite difference solver: Wave and Quake Laboratory—WaveQLab3D, for dynamic earthquake rupture simulations and seismic wave propagation in geometrically complex media.

Traditional SBP finite difference operators used in WaveQLab3D have standard centered finite difference stencils in the interior with one-sided skewed finite stencils close to the boundaries. It is however well known that numerical dispersion properties of centered finite difference operators are less optimal. Therefore, traditional SBP finite difference operators will incur large numerical error for longtime numerical simulations, unless the mesh resolution is increased further.

The goal of this project is to improve the accuracy of WaveQLab3D by replacing the traditional SBP finite difference operators with newly derived and more accurate upwind SBP operators. And benchmark the code against community developed benchmark problems.

Further comments: Programing in Python/Matlab and Fortran 2003 will be required. No previous knowledge of SBP operators is required. This work includes running simulations with WaveQLab3D on a supercomputer. No prior knowledge of that is required.