1 Numerical determination of eigenvalues and eigenvectors

1.1 Preliminaries

Before discussing numerical methods of calculating eigenvalues and eigenvectors we remind you of three results for a matrix $A$ with an eigenvalue $\lambda$ and associated eigenvector $X$ .

• $A^{-1}$ (if it exists) has an eigenvalue $\frac{1}{\lambda }$ with associated eigenvector $X$ .
• The matrix $\left(A-kI\right)$ has an eigenvalue $\left(\lambda -k\right)$ and associated eigenvector $X$ .
• The matrix ${\left(A-kI\right)}^{-1}$ , i.e. the inverse (if it exists) of the matrix $\left(A-kI\right)$ , has eigenvalue $\frac{1}{\lambda -k}$ and corresponding eigenvector $X$ .

  Here $k$ is any real number.

Task!

The matrix $A=\left[\begin{array}{ccc}\hfill 2\hfill & \hfill 1\hfill & \hfill 1\hfill \\ \hfill 1\hfill & \hfill 2\hfill & \hfill 1\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill 5\hfill \end{array}\right]$ has eigenvalues $\lambda =5,3,1$ with associated eigenvectors $\left[\begin{array}{c}\hfill 1∕2\hfill \\ \hfill 1∕2\hfill \\ \hfill 1\hfill \end{array}\right],\left[\begin{array}{c}\hfill 1\hfill \\ \hfill 1\hfill \\ \hfill 0\hfill \end{array}\right],\left[\begin{array}{c}\hfill 1\\ \hfill -1\\ \hfill 0\end{array}\right]$ respectively.

The inverse $A^{-1}$ exists and is $A^{-1}=\frac{1}{3}\left[\begin{array}{ccc}\hfill 2& \hfill -1& \hfill -5\\ \hfill -1& \hfill 2& \hfill -5\\ \hfill 0& \hfill 0& \hfill \frac{3}{5}\end{array}\right]$

Without further calculation write down the eigenvalues and eigenvectors of the following matrices:

1. $A^{-1}$
2. $\left[\begin{array}{ccc}\hfill 3\hfill & \hfill 1\hfill & \hfill 1\hfill \\ \hfill 1\hfill & \hfill 3\hfill & \hfill 1\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill 6\hfill \end{array}\right]$
3. ${\left[\begin{array}{ccc}\hfill 0\hfill & \hfill 1\hfill & \hfill 1\hfill \\ \hfill 1\hfill & \hfill 0\hfill & \hfill 1\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill 3\hfill \end{array}\right]}^{-1}$

Answer

1.2 The power method

This is a direct iteration method for obtaining the dominant eigenvalue (i.e. the largest in magnitude), say ${\lambda }_1$ , for a given matrix $A$ and also the corresponding eigenvector.

We will not discuss the theory behind the method but will demonstrate it in action and, equally importantly, point out circumstances when it fails .

Task!

Let $A=\left[\begin{array}{cc}\hfill 4\hfill & \hfill 2\hfill \\ \hfill 5\hfill & \hfill 7\hfill \end{array}\right]$ . By solving $det\; <!--nolimits\; -->\left(A-\lambda I\right)=0$ obtain the eigenvalues of $A$ and also obtain the eigenvector associated with the dominant eigenvalue.

Answer

We shall now demonstrate how the power method can be used to obtain ${\lambda }_1=9$ and $X=\left[\begin{array}{c}\hfill 0.4\hfill \\ \hfill 1\hfill \end{array}\right]$ where $A=\left[\begin{array}{cc}\hfill 4\hfill & \hfill 2\hfill \\ \hfill 5\hfill & \hfill 7\hfill \end{array}\right]$ .

• We choose an arbitrary $2\times 1$ column vector

  $X^{\left(0\right)}=\left[\begin{array}{c}\hfill 1\hfill \\ \hfill 1\hfill \end{array}\right]$

• We premultiply this by $A$ to give a new column vector $X^{\left(1\right)}$ :

  $X^{\left(1\right)}=\left[\begin{array}{cc}\hfill 4\hfill & \hfill 2\hfill \\ \hfill 5\hfill & \hfill 7\hfill \end{array}\right]\left[\begin{array}{c}\hfill 1\hfill \\ \hfill 1\hfill \end{array}\right]=\left[\begin{array}{c}\hfill 6\\ \hfill 12\end{array}\right]$

• We ‘normalize’ $X^{\left(1\right)}$ to obtain a column vector $Y^{\left(1\right)}$ with largest component 1: thus

  $Y^{\left(1\right)}=\frac{1}{12}\left[\begin{array}{c}\hfill 6\\ \hfill 12\end{array}\right]=\left[\begin{array}{c}\hfill 1∕2\hfill \\ \hfill 1\hfill \end{array}\right]$

• We continue the process

  $X^{\left(2\right)}=A{Y}^{\left(1\right)}=\left[\begin{array}{cc}\hfill 4\hfill & \hfill 2\hfill \\ \hfill 6\hfill & \hfill 7\hfill \end{array}\right]\left[\begin{array}{c}\hfill 1∕2\hfill \\ \hfill 1\hfill \end{array}\right]=\left[\begin{array}{c}\hfill 4\hfill \\ \hfill 9.5\hfill \end{array}\right]$

  $Y^{\left(2\right)}=\frac{1}{9.5}\left[\begin{array}{c}\hfill 4\hfill \\ \hfill 9.5\hfill \end{array}\right]=\left[\begin{array}{c}\hfill 0.421053\hfill \\ \hfill 1\hfill \end{array}\right]$

Task!

Continue this process for a further step and obtain $X^{\left(3\right)}$ and $Y^{\left(3\right)}$ , quoting values to 6 d.p.

Answer

The first 8 steps of the above iterative process are summarized in the following table (the first three rows of which have been obtained above):

Table 1

In Table 1, ${\alpha }_r$ refers to the largest magnitude component of $X^{\left(r\right)}$ which is used to normalize $X^{\left(r\right)}$ to give $Y^{\left(r\right)}$ . We can see that ${\alpha }_r$ is converging to 9 which we know is the dominant eigenvalue ${\lambda }_1$ of $A$ . Also $Y^{\left(r\right)}$ is converging towards the associated eigenvector ${\left[0.4,1\right]}^T$ .

Depending on the accuracy required, we could decide when to stop the iterative process by looking at the difference $\left|{\alpha }_r-{\alpha }_{r-1}\right|$ at each step.

Task!

Using the power method obtain the dominant eigenvalue and associated

eigenvector of

$A=\left[\begin{array}{ccc}\hfill 3& \hfill -1& \hfill 0\\ \hfill -2& \hfill 4& \hfill -3\\ \hfill 0& \hfill -1& \hfill 1\end{array}\right]$ using a starting column vector $X^{\left(0\right)}=\left[\begin{array}{c}\hfill 1\hfill \\ \hfill 1\hfill \\ \hfill 1\hfill \end{array}\right]$

Calculate $X^{\left(1\right)},$ $Y^{\left(1\right)}$ and ${\alpha }_1$ :

Answer

Carry out the next two steps of this iteration to obtains $X^{\left(2\right)},Y^{\left(2\right)},{\alpha }_2$ and $X^{\left(3\right)},Y^{\left(3\right)},{\alpha }_3$ :

Answer

After just 3 iterations there is little sign of convergence of the normalizing factor ${\alpha }_r$ . However the next two values obtained are

${\alpha }_4=5.7347{\alpha }_5=5.4774$

and, after 14 iterations, $\left|{\alpha }_{14}-{\alpha }_{13}\right|<0.0001$ and the power method converges, albeit slowly, to

${\alpha }_{14}=5.4773$

which (correct to 4 d.p.) is the dominant eigenvalue of $A$ . The corresponding eigenvector is

$\left[\begin{array}{c}\hfill -0.4037\hfill \\ \hfill 1\hfill \\ \hfill -0.2233\hfill \end{array}\right]$

It is clear that the power method requires, for its practical execution, a computer.

1.3 Problems with the power method

1. If the initial column vector $X^{\left(0\right)}$ is an eigenvector of $A$ other than that corresponding to the dominant eigenvalue, say ${\lambda }_1$ , then the method will fail since the iteration will converge to the wrong eigenvalue, say ${\lambda }_2$ , after only one iteration (because $A{X}^{\left(0\right)}={\lambda }_2{X}^{\left(0\right)}$ in this case).
2. It is possible to show that the speed of convergence of the power method depends on the ratio

   $\frac{\text{ magnitude of dominant eigenvalue}{\lambda }_1}{\text{ magnitude of next largest eigenvalue}}$

   If this ratio is small the method is slow to converge.

   In particular, if the dominant eigenvalue ${\lambda }_1$ is complex the method will fail completely to converge because the complex conjugate ${\bar{\lambda }}_1$ will also be an eigenvalue and $\left|{\lambda }_1\right|=\left|{\bar{\lambda }}_1\right|$

3. The power method only gives one eigenvalue, the dominant one ${\lambda }_1$ (although this is often the most important in applications).

1.4 Advantages of the power method

1. It is simple and easy to implement.
2. It gives the eigenvector corresponding to ${\lambda }_1$ as well as ${\lambda }_1$ itself. (Other numerical methods require separate calculation to obtain the eigenvector.)

1.5 Finding eigenvalues other than the dominant

We discuss this topic only briefly.

1.5.1 Obtaining the smallest magnitude eigenvalue

If $A$ has dominant eigenvalue ${\lambda }_1$ then its inverse $A^{-1}$ has an eigenvalue $\frac{1}{{\lambda }_1}$ (as we discussed at the beginning of this Section.) Clearly $\frac{1}{{\lambda }_1}$ will be the smallest magnitude eigenvalue of $A^{-1}$ . Conversely if we obtain the largest magnitude eigenvalue, say ${\lambda }_1^{\prime }$ , of $A^{-1}$ by the power method then the smallest eigenvalue of $A$ is the reciprocal, $\frac{1}{{\lambda }_1^{\prime }}.$

This technique is called the inverse power method .

Example

If $A=\left[\begin{array}{ccc}\hfill 3& \hfill -1& \hfill 0\\ \hfill -2& \hfill 4& \hfill -3\\ \hfill 0& \hfill -1& \hfill 1\end{array}\right]$ then the inverse is $A^{-1}=\left[\begin{array}{ccc}\hfill 1\hfill & \hfill 1\hfill & \hfill 3\\ \hfill 2\hfill & \hfill 3\hfill & \hfill 9\\ \hfill 2\hfill & \hfill 3\hfill & \hfill 10\end{array}\right]$ .

Using $X^{\left(0\right)}=\left[\begin{array}{c}\hfill 1\hfill \\ \hfill 1\hfill \\ \hfill 1\hfill \end{array}\right]$ in the power method applied to $A^{-1}$ gives ${\lambda }_1^{\prime }=13.4090$ . Hence the smallest magnitude eigenvalue of $A$ is $\frac{1}{13.4090}=0.0746$ . The corresponding eigenvector is $\left[\begin{array}{c}\hfill 0.3163\hfill \\ \hfill 0.9254\hfill \\ \hfill 1\hfill \end{array}\right].$

In practice, finding the inverse of a large order matrix $A$ can be expensive in computational effort. Hence the inverse power method is implemented without actually obtaining $A^{-1}$ as follows.

As we have seen, the power method applied to $A$ utilizes the scheme:

$X^{\left(r\right)}=A{Y}^{\left(r-1\right)}r=1,2,\dots$

where $Y^{\left(r-1\right)}=\frac{1}{{\alpha }_{r-1}}{X}^{\left(r-1\right)}$ ,   ${\alpha }_{r-1}$ being the largest magnitude component of $X^{\left(r-1\right)}$ .

For the inverse power method we have

$X^{\left(r\right)}=A^{-1}{Y}^{\left(r-1\right)}$

which can be re-written as

$A{X}^{\left(r\right)}=Y^{\left(r-1\right)}$

Thus $X^{\left(r\right)}$ can actually be obtained by solving this system of linear equations without needing to obtain $A^{-1}$ . This is usually done by a technique called $LU$ decomposition i.e. writing $A$ (once and for all) in the form

$A=LU$ $L$ being a lower triangular matrix and $U$ upper triangular.

1.5.2 Obtaining the eigenvalue closest to a given number $p$

Suppose ${\lambda }_k$ is the (unknown) eigenvalue of $A$ closest to $p$ . We know that if ${\lambda }_1,{\lambda }_2,\dots ,{\lambda }_n$ are the eigenvalues of $A$ then ${\lambda }_1-p,{\lambda }_2-p,\dots ,{\lambda }_n-p$ are the eigenvalues of the matrix $A-pI$ . Then ${\lambda }_k-p$ will be the smallest magnitude eigenvalue of $A-pI$ but $\frac{1}{{\lambda }_k-p}$ will be the largest magnitude eigenvalue of ${\left(A-pI\right)}^{-1}$ . Hence if we apply the power method to ${\left(A-pI\right)}^{-1}$ we can obtain ${\lambda }_k$ . The method is called the shifted inverse power method.

1.5.3 Obtaining all the eigenvalues of a large order matrix

In this case neither solving the characteristic equation $det\; <!--nolimits\; -->\left(A-\lambda I\right)=0$ nor the power method (and its variants) is efficient.

The commonest method utilized is called the QR technique . This technique is based on similarity transformations i.e. transformations of the form

$B=M^{-1}AM$

where $B$ has the same eigenvalues as $A$ . (We have seen earlier in this Workbook that one type of similarity transformation is $D=P^{-1}AP$ where $P$ is formed from the eigenvectors of $A$ . However, we are now, of course, dealing with the situation where we are trying to find the eigenvalues and eigenvectors of $A$ .)

In the $QR$ method $A$ is reduced to upper (or lower) triangular form. We have already seen that a triangular matrix has its eigenvalues on the diagonal.

For details of the $QR$ method, or more efficient techniques, one of which is based on what is called a Householder transformation, the reader should consult a text on numerical methods.