Algorithm 840: computation of grid points, quadrature weights and derivatives for spectral element methods using prolate spheroidal wave functions---prolate elements

• Authors: John P. Boyd
• Venue: ACM Transactions on Mathematical Software (TOMS)
• Year: 2005

High order domain decomposition methods using a basis of Legendre polynomials, known variously as “spectral elements” or “p-type finite elements,” have become very popular. Recent studies suggest that accuracy and efficiency can be improved by replacing Legendre polynomials by prolate spheroidal wave functions of zeroth order. In this article, we explain the practicalities of computing all the numbers needed to switch bases: the grid points xj, the quadrature weights wj, and the values of the prolate functions and their derivatives at the grid points. The prolate functions themselves are computed by a Legendre-Galerkin discretization of the prolate differential equation; this yields a symmetric tridiagonal matrix. The prolate functions are then defined by Legendre series whose coefficients are the eigenfunctions of the matrix eigenproblem. The grid points and weights are found simultaneously through a Newton iteration. For large N and c, the iteration diverges from a first guess of the Legendre-Lobatto points and weights. Fortunately, the variations of the xj and wj with c are well-approximated by a symmetric parabola over the whole range of interest. This makes it possible to bypass the continuation procedures of earlier authors.