Tom Hagstrom's Radiation Boundary Condition Page

This purpose of this page is to collect and disseminate useful information related to the accurate and efficient near-field truncation of the computational domain for the simulation of wave propagation problems in the time domain. A multitude of methods have been proposed for this task; see for example the relatively recent review articles:

  • New results on absorbing layers and radiation boundary conditions , Topics in Computational Wave Propagation, M. Ainsworth et al eds., Springer-Verlag, 2003, 1-42.
  • Radiation boundary conditions for Maxwell's equations: a review of accurate time-domain formulations , J. Comput. Math., 25, 2007, 305-336.
    Or the recent lecture:
  • A Lecture at the GAMM Workshop 2009

    Current contents are limited to two techniques which can be used for problems such as the scalar wave equation, Maxwell's equations and the acoustic system in uniform, isotropic media. The first involves compressions of the radiation boundary kernel arising in nonlocal formulations of exact conditions on special boundaries. The second involves optimal local boundary condition sequences applicable on polygonal boundaries.


    Exact formulas for radiation boundary conditions on special boundaries - planes, cylinders, and spheres - can be obtained by separation of variables. For the standard models these formulas all involve a nonlocal operator of the form F-1KF where F is a spatial harmonic transform (Fourier on the plane, spherical harmonics on the sphere) and K is a temporal convolution. The direct evaluation of K is inefficient as it involves the entire time history of the solution on the boundary. However, a fast, low-memory approximate evaluation is achieved by replacing the convolution kernel by a sum of exponentials and noting that convolution with an exponential kernel is equivalent to solving an ode and thus requires no global memory. It can be proven that the kernels arising from the scalar wave equation and related models can be approximated to an accuracy ε using O(log 1/ε· log cT/λ) exponentials where T is the simulation time and λ is the wavelength. See:

  • Rapid evaluation of nonreflecting boundary kernels for time-domain wave propagation , SIAM J. Numer. Anal., 37, 2000, 1138-1164.
  • Nonreflecting boundary conditions for the time-dependent wave equation , J. Comput. Phys., 180, 2002, 270-296.

    To use these compressed approximations one only needs to know the amplitudes and exponents in the sum-of-exponential approximations. You can access these here for ε=1E-6:

  • Description of the data files
  • Cylinder kernel
  • Plane kernel
  • Sphere kernel


    There are two drawbacks to the nonlocal conditions. The first is the need to use spatial harmonic transforms, which is some expense in 3+1 dimensions and involves some effort to couple with the interior scheme. The second, and most important, is the restriction on the shape of the artificial boundary. Particularly for high-aspect-ratio scatterers it would be more efficient to bound the computational domain by a box rather than a sphere. This can be done using local methods such as local radiation boundary condition sequences or perfectly matched absorbing layers. However, these local methods can suffer severe accuracy degradation over time. To correct this defect we have introduced a new parametrization of local boundary condition sequences which we call Complete Radiation Boundary Conditions (CRBCs). These involve involve inhomogeneous rational approximants to the transform of the exact planar kernel which interpolate it in the right half-plane. Boundary conditions are built out of modified ``Higdon-type'' operators: ajd/dt+(1-aj2)/(aj T) ±cd/dn. The cosines aj are determined by specifying the boundary condition order and the dimensionless parameter η=δ/cT where T is the simulation time, c is the wavespeed, and δ is the minimal separation between the artificial boundary and any sources, scatterers, or other inhomogeneities. For details see the preprint:

  • Complete radiation boundary conditions: minimizing the long time error growth of local methods.

    A table of parameters aj for boundary condition orders P=4,8,... (tolerances > 1E-8) and η=1E-2,...,1E-6 along with maximum values of the reflection coefficient may be accessed below. Note that the number of cosines required by a condition of order P is 2P+2.

  • Optimal cosines for certain values of P and η.

    These cosines are computed by the following MATLAB function which implements the Remez algorithm. Its inputs are η and P and the outputs are the cosines (aj, j=1, ... , 2P+2) and the maximum of the complex reflection coefficient. Note that the function may fail due to conditioning issues if P is chosen too large.

  • A MATLAB function for computing optimal cosines.

    We note that direct applications to second-order formulations have so far used different parametrizations based on a combination of Gauss-Lobatto and Yarvin-Rokhlin quadrature nodes. For details see:

  • Radiation boundary conditions for time-dependent waves based on complete plane wave expansions, J. Comput. Appl. Math., to appear.


    I greatfully acknowledge the contribution of many collaborators in this effort including Brad Alpert, Dan Givoli, Leslie Greengard and Tim Warburton. The work is currently suppoorted by the National Science Foundation via grant DMS-06010067 and has also been supported in part by ARO Grant DAAD19-03-1-0146, AFOSR Contract FA9550-05-1-0473, and the Israel-US Binational Science Foundation. Any conlusions or recommendations expressed here are my own and do not necessarily reflect the views of NSF, ARO, AFOSR, BSF, or my collaborators. If you have any questions, comments, or suggestions please contact me at:
    thagstrom at smu dot edu
    Last updated: January 5, 2009.