Finding the Eigenvalues of a Matrix

Computes extremal eigenvalues of a matrix

Description

The program shows how one can find extremal eigenvalues (the largest and the smallest) as well as the eigenvalue nearest to some target value. These methods are described in great details in the book by Kenneth J. Beers, Numerical Methods for Chemical Engineering, Applications in Matlab, Cambridge University Press, 2007.

• MATLAB release: MATLAB 7.0.1 (R14SP1)

[Scott Miller]

Neat. Most people will want to package this in a function with error checking, etc., but this is still a nice technique (especially the part about "nearest" eigenvalue/vector). The "neat" factor for me outweighs the lack of modularity by far.

[John D'Errico]

This code shows how to use the power method for eigenvalues. While I have no problem with a code that teaches the power method (instead of just using eig instead) it should also explain where the power method will fail.

For example, this code will have problems if the target eigenvalue is complex, since the power method will always predict a real eigenvalue for a real matrix. (Note that many matrices in engineering applications will be symmetric, and even positive definite, but note that no such requirement is stated.)

Eigenvalue pairs that are close in magnitude will also cause problems - a small fixed number of iterations will be completely inadequate there.

Next, this method may fail completely for some problems. Consider this trivial example:

A = diag([3 4 5])
A =
     3 0 0
     0 4 0
     0 0 5

Yes, we know the eigenvalues quite well, so its a good test case. 4 is indeed exactly an eigenvalue. Finding the eigenvalue closest to 4 does not work as well as you might wish:

v=[1.,1,1]';
for i=2:10
    v=inv(A-4*eye(3))*v/norm(inv(A-4*eye(3))*v);
end
lamda3=1/(v'*inv(A-4*eye(3))*v)+4

Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.
Warning: Matrix is singular to working precision.

lamda3 =
   NaN

Next, what happens when there are two eigenvalues with the same magnitude? Find the "largest" one:

A = diag([-1 1 .5]);

v=[1.,1,1]';
for i=2:100
    v=A*v/norm(A*v);
end
lamda1=v'*A*v

lamda1 =
   2.2371e-17

Nope. That was not it. Not even close.

Also beware of using a vector of ones as the starting eigenvector. Find the largest eigenvalue of this matrix using the methods in this submission:
 
A = [3 -1 -2;-1 3 -2;-2 -2 4];

v=[1.,1,1]';
for i=2:10
    v=A*v/norm(A*v);
end
lamda1=v'*A*v

Warning: Divide by zero.

lamda1 =
   NaN

The eigenvalues are actually 0, 4, and 6. Lets try a different starting vector:

v=[1, -1, 0]';
for i=2:10
    v=A*v/norm(A*v);
end
lamda1=v'*A*v

lamda1 =
            4

Nope. It still failed. Eig works nicely though.

eig(A)
ans =
  -4.4409e-16
            4
            6

The point is, there are many reasons to avoid the power method. If this code is to be of any value, it should suggest when to be careful.

[Gurpreet Singh]

This method works for only 3x3 matrix. Can you generate a program for nxn matrix that means any size?