David's Blog

Living a quiet life in Coquitlam, B.C.

Name:
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.

References:

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.
      1976

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

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: , , ,

Sunday, July 18, 2010

Eigenvalue and Eigenvector Calculator

The web page to which the title of this post is hot linked is the most recent utility I have added for computing the eigenvalues and eigenvectors of a real general matrix. It accepts square matrices of up to size twelve by twelve (this limit was arbitrarily picked by me, because that was the largest size matrices with which I worked during my undergraduate engineering studies).

As cited on the page itself, this eigenvalue and eigenvector solver is a translation of routines that come from the EISPACK collection of routines, specifically, balanc, hqr2, elmhes, eltran, and balbak.

Note that the results of the computation may be complex, but the inputs are all assumed to be real entries.

Also note that this utility can solve symmetric matrices, but using a program specifically written for symmetric matrices would be more efficient.

Labels: , , , , , , ,

Sunday, June 06, 2010

Online Matrix Calculator

The web page to which the title of this post is hot-linked is the most recent utility I have added for solving simultaneous equations. It features a few changes and additions compared to the simultaneous equation solvers already posted on the site.

First of all, as mentioned in a previous post, the method of accepting input has been changed; instead of accepting input from an individual text box for each entry, data is now input all at once via a textarea box. Similarly, instead of outputting results to individual text boxes, results are output to the body of the web page. These changes are meant to make it more convenient to copy data to and from the utility, thereby reducing the possibility of introducing errors into the input and output. Note that the output can look confusing for larger matrices; lines wrap to stay within the body of the web page, breaking up rows of the matrix. For example, in the case of a 12 x 12 matrix, on my computer each row of the matrix is output on three lines. However, if the entire output section is copied to a plaintext document, the line wraps go back to where they belong (i.e. - the matrix displays as expected: twelve lines of twelve entries each).

The program, itself, offers some additional capabilities compared to the other simultaneous equation solvers on this site: in addition to being able to solve larger systems (up to twelve equations in twelve unknowns), this utility computes the determinant, the LU-Decomposition, and the inverse of the [A] matrix. In fact, much of this additional capability is included in the existing programs, but not used. The algorithm underlying the program for solving N Equations in N Unknowns performs an LU-Decomposition as one of the steps toward solving the system, but the results of the LU-Decomposition are not output by the other Javascript pages. Furthermore, once an LU-Decomposition has been performed on a matrix, computing the determinant and inverse of that matrix is quite simple. These matrix properties, too, are not output by the other Javascript pages presently on site. This latest utility takes advantage of the fact that, since most of the work is already done (the LU-Decomposition), computing and outputting a few additional properties of the [A] matrix is quick and relatively easy. So, this Javascript utility does compute and output those properties:
1) the [P], [L], and [U] matrices are output, where [P] is the permutation matrix (which keeps track of the permutations that were involved in the decomposition), [L] is the Lower Matrix, and [U] is the Upper Matrix;
2) the determinant of the [A] matrix is computed and output; and
3) the inverse of the [A] matrix is computed and output, [A]-1.

Note that this utility may be used to compute these properties of the [A] matrix without solving a valid system of simultaneous equations. For example, say a user only wants to compute the matrix determinant. In this case, the user would enter a dummy b vector. I use (1 1 1 . . . )T, but it doesn’t matter what the b vector is, as long as it has the appropriate number of entries. The utility would then be executed and the results for the x vector ignored. The properties of [A] are not effected by b, so entering a fake b vector doesn’t change the determinant, LU-Decomposition, or inverse of [A].

Labels: , , , , , , , , , , ,