Numerical experiments on the accuracy of the Chebyshev–Frobenius companion matrix method for finding the zeros of a truncated series of Chebyshev polynomials

For a function f(x)f(x) that is smooth on the interval x∈[a,b]x∈[a,b] but otherwise arbitrary, the real-valued roots on the interval can always be found by the following two-part procedure. First, expand f(x)f(x) as a Chebyshev polynomial series on the interval and truncate for sufficiently large N. Second, find the zeros of the truncated Chebyshev series. The roots of an arbitrary polynomial of degree N  , when written in the form of a truncated Chebyshev series, are the eigenvalues of an N×NN×N matrix whose elements are simple, explicit functions of the coefficients of the Chebyshev series. This matrix is a generalization of the Frobenius companion matrix. We show by experimenting with random polynomials, Wilkinson's notoriously ill-conditioned polynomial, and polynomials with high-order roots that the Chebyshev companion matrix method is remarkably accurate for finding zeros on the target interval, yielding roots close to full machine precision. We also show that it is easy and cheap to apply Newton's iteration directly to the Chebyshev series so as to refine the roots to full machine precision, using the companion matrix eigenvalues as the starting point. Lastly, we derive a couple of theorems. The first shows that simple roots are stable under small perturbations of magnitude εε to a Chebyshev coefficient: the shift in the root x*x* is bounded by ε/df/dx(x*)+O(ε2)ε/df/dx(x*)+O(ε2) for sufficiently small εε. Second, we show that polynomials with definite parity (only even or only odd powers of x) can be solved by a companion matrix whose size is one less than the number of nonzero coefficients, a vast cost-saving.