A fast FFT-based discrete Legendre transform

An ${\mathcal O}(N(\log N)^2/\log \log N)$ algorithm for the computation of the discrete Legendre transform and its inverse is described. The algorithm combines a recently developed fast transform for converting between Legendre and Chebyshev coefficients with a Taylor series expansion for Chebyshev polynomials about equally spaced points in the frequency domain. Both components are based on the fast Fourier transform, and as an intermediate step we obtain an ${\mathcal O}(N\log N)$ algorithm for evaluating a degree-$(N-1)$ Chebyshev expansion on an $N$-point Legendre grid. Numerical results are given to demonstrate performance and accuracy.