The orthogonal transformations P_{pq} annihilate element (p,q) of an object matrix.
Successive transformations undo previously set zeros, but the off-diagonal terms eventually get smaller and smaller, until they get nearly zero, and the only thing that's left is the diagonal with eigenvalues.
Taking the product of all the transformations P_{pq} yields the matrix of eigenvectors.
The method is absolutely foolproof for symmetric real matrices. But it is slow. Painfully slow for matrices larger than .
It is a simple and stable algorithm though, and it parallelises well too.
The Jacobi rotation matrix has the form:
(3.5) |
c | = | (3.6) | |
s | = | (3.7) |
c^{2} + s^{2} = 1 | (3.8) |
It is quite easy to see what the effect of equation (3.9) is going to be on selected terms.
First we need to come up with an expression that describes
a generic term of matrix
P_{pq} in terms of
Kronecker deltas:
P_{ij} | = | ||
= | (3.10) |
Now we can evaluate a'_{rp} for, say,
and (assume summation over dummy indexes i and j)
a'_{rp} | = | P^{-1}_{ri} a_{ij} P_{jp} = P_{ir} a_{ij} P_{jp} | |
= | |||
= | |||
= | |||
= | |||
= | (3.11) |
In summary:
a'_{rp} | = | (3.12) | |
a'_{rq} | = | (3.13) | |
a'_{pp} | = | c^{2} a_{pp} + s^{2} a_{qq} - 2 s c a_{pq} | (3.14) |
a'_{qq} | = | s^{2} a_{pp} + c^{2} a_{qq} + 2 s c a_{pq} | (3.15) |
a'_{pq} | = | (3.16) |
The purpose of the Jacobi rotation
P_{pq} is to
kill a'_{pq}. Thus:
(3.17) |
Now, observe that:
(3.19) |
Next recall the following simple algebraic identity:
(3.21) |
Solving equation (3.22) with respect to t yields the following:
Now, this solution can be rewritten in a computationally more convenient
form. Consider the + case:
(3.23) |
Once we have the t, we get c and s as follows^{3.1}
c | = | (3.24) | |
s | = | c t | (3.25) |
Equations (3.12) through (3.16) are now rewritten to
minimize the round off error and to make them look like
the new quantity is equal to the old one plus a small correction.
And so,
(3.27) |
a_{pp}' | = | c^{2} a_{pp} + s^{2} a_{qq} - 2sca_{pq} | |
= | |||
= | |||
= | |||
= | a_{pp} - t a_{pq} | (3.28) |
Similarly from the same equation (3.26) we get:
(3.29) |
a_{qq}' | = | s^{2} a_{pp} + c^{2} a_{qq} + 2 s c a_{pq} | |
= | |||
= | |||
= | |||
= | a_{qq} + t a_{pq} | (3.30) |
For a'_{rp} and a'_{rq} the computation is more trivial:
a'_{rp} | = | c a_{rp} - s a_{rq} | |
= | a_{rp} + (c - 1) a_{rp} - s a_{rq} | ||
= | (3.31) |
a'_{rq} | = | c a_{rq} + s a_{rp} | |
= | a_{rq} + (c - 1) a_{rq} + s a_{rp} | ||
= | (3.32) |
(3.33) |