  We describe a fast, simple, and stable transform of Chebyshev expansion
coefficients to Jacobi expansion coefficients and its inverse based on the
numerical evaluation of Jacobi expansions at the Chebyshev--Lobatto points.
This is achieved via a decomposition of Hahn's interior asymptotic formula into
a small sum of diagonally scaled discrete sine and cosine transforms and the
use of stable recurrence relations. It is known that the Clenshaw--Smith
algorithm is not uniformly stable on the entire interval of orthogonality.
Therefore, Reinsch's modification is extended for Jacobi polynomials and
employed near the endpoints to improve numerical stability.

Numerical Analysis