We can now describe the doubling generalization of Schur's algorithm, which forms the basis for the first phase of our superfast Toeplitz solver.
Recall that if we perform m steps of Schur's algorithm on the Schur function , we can obtain the mth Schur polynomials and . Moreover, these polynomials are determined by the first m coefficients of and . Let and denote the polynomials of degree less than m formed from these coefficients. To describe the doubling procedure, we assume that and have been computed from and , and we seek to compute and from and .
Having obtained and , the mth Schur function is given by
It can be shown  that
so that the first m terms of the numerator and denominator of can be taken from
We can now perform the doubling step: Since is also a Schur function, we can use the same procedure that computed and from and to compute the Schur polynomials and that result from m steps of Schur's algorithm applied to and .
Let denote the LFT that results from m steps of Schur's algorithm applied to . Then we have , so that . Writing this composition in terms of Schur polynomials, we obtain
This discussion is summarized by the following recursive description of the doubling procedure.
The algorithm can be started by performing Step 0 directly by, for example, setting and . More generally, we can use Schur's algorithm to generate for a small value of , and obtain , from the Schur parameters and the recursions (6.1).
Thus, the generalized Schur algorithm consists of various polynomial multiplications and additions performed in a recursive manner. The multiplication of polynomials can be efficiently performed using standard FFT techniques. This results in an algorithm for computing and , where . The Toeplitz system of equations can then be solved by forming from , by Proposition 6.1, and using a Toeplitz inversion formula to perform the solution phase of the algorithm in arithmetic operations.
A detailed analysis of the implementation for Hermitian Toeplitz matrices is given in . Further refinements for real Toeplitz matrices are described in , where it is shown in that the first phase of the superfast Toeplitz solver can be implemented using fewer than real arithmetic operations for a real Toeplitz matrix , where . Since the Levinson-Durbin algorithm requires more than operations, the superfast algorithm has a smaller operation count for . Experimental results presented in  confirm that the two procedures require approximately the same execution time for n=256, and that the relative efficiency of the superfast algorithm quickly increases for larger values of . Moreover, experimental results indicate that there is little or no degradation in the accuracy of the superfast algorithm compared to the Levinson-Durbin algorithm [4,7].
While the experimental results indicate that this superfast algorithm may be as reliable as the Levinson-Durbin algorithm for the first phase of a Toeplitz solver, no stability analysis along the lines of  or  has yet been performed. However, since the Schur polynomials are generalizations of the Schur parameters, the superfast algorithm is closely related to the hybrid algorithm mentioned at the end of Section 4. This is encouraging in that the superfast algorithm may behave more like the hybrid algorithm than the Levinson-Durbin algorithm.