A stabilized superfast solver for nonsymmetric Toeplitz systems

We present a stabilized superfast solver for nonsymmetric Toeplitz systems Tx=b. An explicit formula for T-1 is expressed in such a way that the matrix-vector product T^-1b can be calculated via FFTs and Hadamard products. This inversion formula involves certain polynomials that can be computed by solving two linearized rational interpolation problems on the unit circle. The heart of our Toeplitz solver is a superfast algorithm to solve these interpolation problems. To stabilize the algorithm, i.e., to improve the accuracy, several techniques are used: pivoting, iterative improvement, downdating, and giving "difficult" interpolation points an adequate treatment. We have implemented our algorithm in Fortran 90. Numerical examples illustrate the effectiveness of our approach.