Pseudospectral time-domain method

2D PSTD simulation of a ground-penetrating radar.

In computational physics, pseudospectral time-domain method (PSTD) is a numerical analysis technique for simulating wave propagation. Being a part of the general class of pseudo-spectral methods, it is an extension of the finite-difference time-domain method (FDTD): in PSTD, the spatial derivative terms of the wave equation are evaluated in the spectral domain using orthogonal bases, such as Fourier or Chebyshev bases while the problem is temporally discretized.

Introduced in the late 1990s, PSTD is widely used in acoustic, electromagnetic and geophysical simulations of complex and large-scale media, with a wide variety of applications such as photoacoustic imaging, geophysical imaging, light scattering, plasma physics and ultrasonics. The method is implemented in open-source acoustics codes such as k-Wave and openPSTD.

Theory

The main concepts behind PSTD can be illustrated through the one-dimensional scalar wave equation:

$\frac{{\partial }^{2}u(x,t)}{\partial x^2}+\frac{1}{v^2}\frac{{\partial }^{2}u(x,t)}{\partial t^2}=0$

where $u(x,t)$ is the field component and $v$ is the wave speed. At a fixed point in time, this field component can be decomposed into an orthogonal basis, such as the Fourier basis. The main concept behind the PSTD algorithm is evaluating the spatial derivatives in this spectral-domain, while using finite difference approximations to compute temporal derivatives. As an example, the derivatives can be represented in the Fourier basis as:

$\frac{\partial u(x,t)}{\partial x}=\frac{\partial }{\partial x}\sum _{n=-\mathrm{\infty }}^{\mathrm{\infty }}{\tilde{u}}_n(t)e^{i{k}_{n}x}=\sum _{n=-\mathrm{\infty }}^{\mathrm{\infty }}i{k}_n{\tilde{u}}_n(t)e^{i{k}_{n}x}={ℱ}^{-1}\left[ikℱ[u(x,t)]\right]$

where $ℱ$ is the Fourier operator and $k$ is the spatial frequency. Fourier transforms can be calculated efficiently through the fast Fourier transform; this formulation enables the spatial resolution of field components in coarser meshes compared to the FDTD method, being only limited by the Nyquist–Shannon sampling theorem. This imposes periodic boundary conditions at the simulation boundaries, effectively causing the propagating fields to wrap around. These effects can be alleviated through the use of perfectly matched layers for the absorption of incident waves. One of the main disadvantages of this approach is the possible artifacts, due to Gibbs phenomenon, which are observed for spatial discontinuities and large contrasts, as well as for discontinuous source conditions. The resolution in this case can be improved with non-uniform fast Fourier transform. Another alternative approach is using Chebyshev polynomials as spatial basis functions, which are evaluated at Chebyshev nodes to minimize Runge phenomenon. This method is generally known as Chebyshev PSTD. Chebyshev PSTD methods can be extended to multidomain decomposition for modeling wave propagation in curved geometries.

Time integration im PSTD algorithms are typically performed through the leapfrog integration scheme as in FDTD; nevertheless, Courant number of PSTD method is $\pi$ times smaller than FDTD, requiring finer time steps per unit mesh to be numerically stable. As an alternative, Runge–Kutta methods can be used to perform time integration.

See also

• Discontinuous Galerkin method
• Finite-difference frequency-domain method
• Spectral element method

Further reading

• Template:Cite book
• Template:Cite book