Talk:Jacobi eigenvalue algorithm

Algorithm Calculation and Example are way off[edit]

The algorithms use of e to supposedly keep track of the eigenvalue approximations is baffling and unneeded. The approximation of the eigenvalues are of course just the diagonal elements of the transformed matrix, which is also the proper values to use for calculating y. Furthermore, the example given is very wrong. The eigenvalues are orders of magnitude too high (as can be verified by simply taking the trace of the original matrix), and the eigenvectors are incorrect. — Preceding unsigned comment added by 88.77.167.5 (talk) 15:29, 31 October 2014 (UTC)[reply]

^ Thanks to whoever posted this note here. I was using the example on this page to test my code and couldn't for the life of me figure out what was going wrong... — Preceding unsigned comment added by 108.203.0.241 (talk) 13:17, 13 January 2015 (UTC)[reply]

capital Gamma[edit]

What does the capital Gamma in the convergence section mean? selfworm^Talk) 00:20, 26 February 2013 (UTC)[reply]

Bug![edit]

Hi all! I've just tried to use the algorithm, but there's a BUG! The maxind method is totaly wrong. A correct one would be:
int maxind(int k)
{
int m = 1;
for (int i = 2; i <= n; i++)
if(abs(S[k][i]) > abs(S[k][m])) m = i;
return m;
}
I hope there are no more errors...

Fixed

Pivoting strategy?[edit]

As far as I know, the Jacobi rotations are not only applied to rows, but also to columns. This means that two entries in each row are changed, and if this happens to eliminate the original maximal elements, we would have to search the entire rows for the next maximal elements. In short: I doubt that the given algorithm always finds the true maximal elements.
The only solution I can think of at the moment would be to put all elements into a heap structure and update the heap property after each elimination step. This would like to a complexity in $O(n^{2}log(n))$ per step. —Preceding unsigned comment added by Sboerm (talk • contribs) 20:51, 15 May 2011 (UTC)[reply]

Quality[edit]

I touched this article some time back, but I did not write it and do not vouch for it. When I discovered it, the article was a complete disaster; I fiddled with the formatting alone, which is not enough. Frankly, it might be better to delete it and start over. (Mea culpa; if I had left it ugly that might have happened already!) --KSmrq^T 13:43, 15 June 2007 (UTC)[reply]

Could someone say (or even better, add to the article), what function $\Gamma (S)$ is in the 'Convergence' section? It's not clear from context, a function that I haven't seen used and hard to search for. This is almost surely the Frobenius norm mentioned earlier in the article. Could someone confirm before I go and change?

Trefethen en Bau, numerical linear algebra mentions that the sum of squares of the off-diagonal entries decreases by at least the factor 1-2/(m^2-m) at each step. (for an mxm matrix) However, it leaves the proof as an exercise to the reader.

I suspect that the Gamma function is supposed to be the sum of squares of the off-diagonal elements. Mastomer (talk) 19:17, 12 August 2013 (UTC)[reply]

Links[edit]

This article should be linked to the "Jacobi Rotation" article.

192.249.47.195 (talk) 14:59, 1 October 2009 (UTC)[reply]

Incorrect calculation of theta[edit]

I'm not at all famailiar with this algorithm, but there appears to be an error in the derivation of theta for a particular i & j.

Given that 0 = cos(2θ)Aij + 0.5sin(2θ)(Aii-Ajj),

tan(2θ) = 2Aij/(Aii - Ajj) is wrong. It should be: tan(2θ) = -2Aij/(Aii - Ajj) or, tan(2θ) = 2Aij/(Ajj - Aii)

I'd like some input before editing the article, however. Snrrub (talk) 08:24, 6 July 2009 (UTC)[reply]

I don't really know a lot about it either, though you seem to be correct. Looking at Gentle's books (such as this one, but he uses exactly the same material in other books), which is the same as the current version of the article. However he also makes other mistakes (claiming that only the 4 elements change), so perhaps we should find a better reference. —3mta3 (talk) 09:22, 6 July 2009 (UTC)[reply]

One of the links in the article[1] gives the same form as you, so we might as well change it.—3mta3 (talk) 09:25, 6 July 2009 (UTC)[reply]

I think it depends on which convention you use for Givens rotations. The article on Givens rotations (which is linked to by this article) has the -s in the upper triangular part, and after some calculations, it seems that in that case, the old version was correct. The article you linked indeed has -s in the lower triangular part, in which case the current version is correct. I definitely think that if we keep the current convention, we should at least mention that we use a different one from the one in the Givens rotation article.81.71.59.155 (talk) 12:39, 25 September 2011 (UTC)[reply]

I can confirm this is an error by comparing with eq 10.24 in Ben Noble first (1969) edition and by trying it with a simple example. Only the version tan(2θ) = -2Sij/(Sii - Sjj) zeros Aij.

Also I think this assumes j>i; i.e. Sij is above the diagonal (and the same θ zeros Sji)

Otherwise you might think that to zero S(j,i) you switch i and j, which would give the negative of the same angle, which zeros neither S(j,i) nor Sij.

If i > j then the formula in the article would be correct. Perhaps also if the transpose of the Givens rotation is used, as mentioned by 81.71.59.155. Eaberry (talk) 17:46, 15 March 2024 (UTC)[reply]

No need for keeping track of changes to determine when to halt algorithm[edit]

In addition to the unnecessary and incorrect use of the e vector to track approximate eigenvalues, the algorithm in the main article unnecessarily keeps track of when an eighenvalue was last changed. All that is necessary to determine when to stop the algorithm is to keep track of the sum of squares of the main diagonal of the transformed matrix. When the difference from one iteration to the next does not change the sum of squares on the main diagonal by very much, the algorithm can halt. Better yet, compare the sum of squares of the diagonal components of the image matrix at each point to the square of the Frobenius norm of the original, this gives an objective score indicating how close to a pure diagonal matrix one has achieved.

Note that this also deals with the infinite-loop case of diagonal matrices, etc. 178.6.80.219 (talk) 11:41, 1 November 2014 (UTC)[reply]

I had the same thought as well, and it can be even simpler when using the classical method of choosing the largest off-diagonal to pivot on (as the example algorithm does) by simply terminating if that largest value is 0. Conceptually, it seems quite silly to keep this extra state when you already have all the information you need once you find the pivot value.

However, while can terminate without the change list, it will typically take longer to do so. This is because of how floating point math works: it can accurately represent values very close to zero, but adding a very small value to a larger value may have no change. For example, 100 + 1e-20 == 100 when using double precision numbers. As a result, while the off-diagonal values will get closer and closer to zero each iteration, they may cease to impact the eigenvalues well before they actually reach zero. This holds true even when using the square of the value for the zero check. (e.g. squaring 1e-20 would yield 1e-40, which is still well within the range for a double precision number close to zero)

When implementing the described algorithm in C with doubles, tracking the changed values converges after 18 iterations. However, when only checking the square of the max pivot value or if some of the intermediate values become sub-normal or zero, it takes 29 iterations to detect that it converges. 174.66.156.190 (talk) 00:49, 30 December 2023 (UTC)[reply]

Algorithm section[edit]

Can someone take the pseudocode in the Algorithm section and translate it into an actual piece of code? It seems to be written in some weird combination of math symbols and several different (old!) languages. Why put the reader through the pain of having to translate it when it could be an actual code snippet instead? Any standard language would be fine. As it is, it's too messy to be easily understood, I can't make heads or tails of it.65.78.19.223 (talk) 01:34, 24 December 2020 (UTC)[reply]

It is called pseudocode because it is not a specific programming language. In my opinion, the pseudocode is easy to read and understand. I am strongly opposed to substituting for a particular language. --DaBler (talk) 09:55, 24 December 2020 (UTC)[reply]