David's Blog

Living a quiet life in Coquitlam, B.C.

Location: Coquitlam, British Columbia, Canada

Sunday, June 10, 2018

Eigenvalues Calculator of a Real Symmetric Matrix

This program computes the eigenvalues of a real symmetric matrix. The program accepts input of a full matrix, A, tests it for symmetry and, if symmetry is confirmed, calculates its eigenvalues. A limit on the size of A is not imposed by the routine. The underlying algorithm is designed to handle real symmetric matrices of any size; however, as the system grows, getting the input data into the DOM textarea on the web page might cause it to hiccup. (Anybody wishing to keep entering larger and larger systems until it breaks, feel free to test it!) The only restriction on A is that its coefficients be real; this program was not written to handle complex coefficients.

The routine is a heavily-edited translation of the EISPACK subprogram RS.F.


Smith, B.T.; J.M. Boyle; J.J. Dongarra; B.S. Garbow; Y. Ikebe; V.C. Klema; and C.B. Moler.
      "Matrix Eigensystem Routines--(EISPACK) Guide"
      Springer-Verlag, Berlin.

Garbow, B.S.; J.M. Boyle; J.J. Dongarra; and C.B. Moler.
      "Matrix Eigensystem Routines--(EISPACK) Guide Extension"
      Springer-Verlag, Berlin.

I would like to emphasize that only eigenvalues are computed; eigenvectors are not. In fact, because only eigenvalues are calculated, only two sub-routines from RS.F were needed: TRED1 and TQL1. The sub-routines for computing both the eigenvalues and eigenvectors (i.e., TRED2 and TQL2) were not needed, so not translated.

TRED1 reduces a real symmetric matrix to a symmetric tridiagonal matrix using orthogonal similarity transformations. It is a translation of the algol procedure tred1, Num. Math. 11, 181-195 (1968) by Martin, Reinsch, and Wilkinson. Handbook for Auto. Comp., vol.ii-linear algebra, 212-226 (1971).

TQL1 finds the eigenvalues of a symmetric tridiagonal matrix by the ql method. It is a translation of the algol procedure tql1, Num. Math. 11, 293-306 (1968) by Bowdler, Martin, Reinsch, and Wilkinson. Handbook for Auto. Comp., vol.ii-linear algebra, 227-240 (1971).

This program is not a straight sub-routine to sub-routine translation: TRED1 and TQL1 are both coded into a single sub-routine and their positions within the program are indicated by comments.

A C++ translation of this program is available on GitHub.
It includes a small sub-routine (Echo2DArray, lines 280 - 291) for outputting the input array, A, to the console. Its call in the main program has been commented out. Anybody wishing to use this program, and to make it smaller, can delete this routine completely.

Note that the main program accepts input of a full symmetric matrix and tests it for symmetry before passing the matrix into the sub-routine that calculates its eigenvalues. Two comments about this aspect of the program:

1) Instead of testing that coefficients across the diagonal are equal, the absolute value of their difference is tested against DBL_EPSILON
(i.e., if (fabs(A[i][j] - A[j][i]) > DBL_EPSILON)).
Testing for true equality (i.e., if (a == b)) when dealing with floating point numbers can lead to unexpected results.

2) The symmetry test is intended to be one more way in which to catch errors: prompting the user to review the data may reveal a typo in a data entry, for instance. The sub-routine that calculates the eigenvalues, in fact, only touches the lower triangular half of the A matrix (coefficients along and below the diagonal.) Coefficients above the diagonal are not touched by the sub-routine.

One final comment before concluding this post. The routine changes the input matrix, A. So, if this routine is incorporated into a larger program, and the original version of A is needed later, be sure to make a copy of A, and then pass the copy into the routine.

Labels: , , ,