Implementing internal interfaces in finite-difference schemes with the Heaviside step function

Abstract

Implementing sharp internal interfaces in finite-difference schemes with high spatial accuracy is challenging. The implementations of interfaces are generally considered accurate to at best second order. The natural way to describe an abrupt change in material parameters is by the use of the Heaviside step function. However, the implementation of the Heaviside step function must be consistent with the discrete sampling on the finite-difference grid. Assuming that the step function takes on the value zero up to some node location and then unity from thereon results in an incorrect wavenumber representation of the Heaviside step function so this representation must be incorrect. However, starting with the proper wavenumber representation of the Heaviside step function and then transforming this spectrum to the space domain give much better accuracy. The interface location appears as a proportionality factor in the phase in the wavenumber domain and can be altered continuously. Thus, the interface can be located anywhere between two node locations. This is a key factor for avoiding stair-case effects from the fields when doing 2D and 3D finite-difference simulations. The proposed method can be used for all systems of partial differential equations that formally can be expressed as a material parameter times a dynamic field on one side of the equal sign and with spatial derivatives on the other side of the equal sign. For geophysical simulations the most important cases will be the Maxwell equations and the acoustic and elastic wave equations.