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 [5] 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 [4]. Further refinements
for real Toeplitz matrices are described in [6],
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
[7] 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 [21] or [11] 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.