Real symmetric tridiagonal eigenproblem: bisection and inverse iteration

There are several ways of finding eigenvalues except for reducing the matrix to tridiagonal form by orthogonal transformations. For example, the number of eigenvalues less than a given number could be easily determined for a symmetric tridiagonal matrix. Let's introduce the function N(w) which is equal to the number of elements of matrix T which are not greater than w. This function is piecewise constant, non-strictly monotone increasing, having discontinuities in points w which are eigenvalues. Thus, the problem of finding eigenvalues is reduced to the problem of finding discontinuities of some piecewise constant function. This can be easily performed by using the bisection of the intervals which contain these discontinuities.

The bisection method allows to find both eigenvalues in a given interval and with given numbers (in ascending order). If it is required to find all eigenvalues, this method is significantly slower than the others, but it can be used to find a small part of eigenvalues.

The bisection method is often used along with the inverse iteration method which allows to find an eigenvector by its corresponding eigenvalue. If it is required to find a small part of a full set of eigenvectors, the inverse iteration is faster than finding all eigenvectors using QL/QR algorithm.

The bisection and inverse iteration methods could be both used on their own and applied to a tridiagonal matrix obtained from a symmetric matrix (this is done by the algorithm which gets a part of eigenvalues and eigenvectors of a symmetric matrix.

Convergence and precision

Let's consider the reliability of such a calculation procedure. This methods are rather reliable when dealing with matrices with well-separated eigenvalues. If some eigenvalues are close to each other, the inverse iteration method can return an inaccurate set of eigenvectors (the closer the eigenvalues are, the less accurate the eigenvectors will be).

It should be understood that such difficulties are typical of all methods of finding eigenvalues and no satisfactory solution has yet been found. Nevertheless, this method is applicable to the majority of real problems. If the algorithm is applied to the problem which it couldn't solve, the subroutine returns False. In this case, it is reasonable to use QL/QR algorithm which is more stable than the inverse iteration method, and then select a required vector from the full set of eigenvectors found (it is rather a slow way but you will get a solution).

Subroutine description

There are two subroutines for finding eigenvalues and eigenvectors. The first, SMatrixTDEVDR, allows to find the eigenvalues from a given interval. The second subroutine, SMatrixTDEVDI, allows to find the eigenvalues (and the corresponding eigenvectors) having given indexes.

This algorithm is transferred from the LAPACK library.

This article is licensed for personal use only.

Download ALGLIB for C++ / C# / ...

ALGLIB Project offers you two editions of ALGLIB:

ALGLIB Free Edition:
delivered for free
offers full set of numerical functionality
extensive algorithmic optimizations
no low level optimizations
non-commercial license

ALGLIB Commercial Edition:
flexible pricing
offers full set of numerical functionality
extensive algorithmic optimizations
high performance (SMP, SIMD)
commercial license with support plan

Links to download sections for Free and Commercial editions can be found below:

ALGLIB 3.12.0 for C++

C++ library.
Delivered with sources.
Monolithic design.
Extreme portability.

ALGLIB 3.12.0 for C#

C# library with native kernels.
Delivered with sources.
VB.NET and IronPython wrappers.
Extreme portability.

ALGLIB 3.12.0 for Delphi

Delphi wrapper around C core.
Delivered as precompiled binary.
Compatible with FreePascal.

ALGLIB 3.12.0 for CPython

CPython wrapper around C core.
Delivered as precompiled binary.