Copyright
©
by SIAM. Unauthorized reproduction of this article is prohibited.
SIAM J. M
ATRIX
A
NAL.
A
PPL
.
c
2010 Society for Industrial and Applied Mathematics
Vol. 31, No. 5, pp. 2531–2552
STABILITY OF THE LEVINSON ALGORITHM FOR
TOEPLITZ-LIKE SYSTEMS
∗
P. FAVATI
†
, G. LOTTI
‡
,
AND
O. MENCHI
§
Abstract.
Numerical stability of the Levinson algorithm, generalized for Toeplitz-like systems,
is studied. Arguments based on the analytic results of an error analysis for ﬂoating point arithmetic
produce an upper bound on the norm of the residual vector, which grows exponentially with respect
to the size of the problem. The base of such an exponential function can be small for diagonally
dominant Toeplitz-like matrices. Numerical experiments show that, for these matrices, Gaussian
elimination by row and the Levinson algorithm have residuals of the same order of magnitude. As
expected, the empirical results point out that the theoretical bound is too pessimistic.
Key words.
Levinson algorithm, Toeplitz-like matrices, stability
AMS subject classifications.
15A06, 15B05
DOI.
10.1137/090753619
1. Introduction.
Toeplitz systems arise frequently in linear algebra (see [2] for
a list of possible sources), and special fast and superfast algorithms have been devised
to solve them. Starting from the original Durbin algorithm to solve the Yule–Walker
equations, the Levinson algorithm has been proposed for symmetric positive defi-
nite Toeplitz matrices and extended to the case of general Toeplitz matrices [8]. The
Levinson algorithm is a fast method, i.e., it has a cost of
O
(
N
2
) operations,
N
being
the size of the system. Unfortunately, when a simple operation like multiplication or
inversion or low rank modification is applied to a Toeplitz matrix, the Toeplitz struc-
ture is lost and more general structures must be considered. The class of Toeplitz-like
matrices, which is closed for the most common operations applied in numerical algo-
rithms, seems ideal from this point of view. It is based on the concept of displacement
rank introduced in [15] and has been studied by many authors (see, for example,
[11, 13, 14]). The displacement operator allows a compact representation of the ma-
trices of this class by means of a set of generators. The standard situation, when
dealing with a structured matrix, assumes that its entries are not exactly known but
can be computed from the generators when they are needed.
For the Toeplitz-like matrices, fast and superfast algorithms have been devised
as well (see the extensive bibliography in [14]), but the question of their stability is
still a matter of discussion. The Levinson algorithm, too, has been generalized for this
class of matrices, maintaining its computational cost [6, 13].
The numerical stability of the original Levinson algorithm for symmetric positive
definite matrices has been proved in [2, 5]. Look-ahead modifications have been pro-
posed in [3, 4] to obviate instability in the general Toeplitz case. We are interested in
studying the stability of the Levinson algorithm generalized for Toeplitz-like systems,