This notebook contains an excerpt from the Python Programming and Numerical Methods - A Guide for Engineers and Scientists, the content is also available at Berkeley Python Numerical Methods.

The copyright of the book belongs to Elsevier. We also have this interactive book online for a better learning experience. The code is released under the MIT license. If you find this content useful, please consider supporting the work on Elsevier or Amazon!

# The QR Method¶

The QR method is a preferred iterative method to find all the eigenvalues of a matrix (but not the eigenvectors at the same time). The idea is based on the following two concepts

1. similar matrices will have the same eigenvalues and associated eigenvectors. Two square matrices $$A$$ and $$B$$ are similar if:

$A = C^{-1}BC$

where $$C$$ is an invertible matrix.

1. The QR method is a way to decompose a matrix into two matrices $$Q$$ and $$R$$, where $$Q$$ is an orthogonal matrix, and $$R$$ is an upper triangular matrix. An orthogonal matrix has the features: $$Q^{-1} = Q^T$$, which means $$Q^{-1}Q=Q^TQ=I$$.

How do we link these two concepts to find the eigenvalues? Say we have a matrix $$A_0$$ whose eigenvalues must be determined. At the $$k$$-th step (starting with $$k = 0$$), we can perform the QR decomposition and get $$A_k = Q_kR_k$$, where $$Q_k$$ is an orthogonal matrix and $$R_k$$ is an upper triangular matrix. We then form $$A_{k+1} = R_kQ_k$$, which we note that

$A_{k+1} = R_kQ_k = Q^{-1}_kQ_kR_kQ_k = Q^{-1}_kA_kQ_k$

therefore, all the $$A_k$$ are similar, as we discussed above they all have the same eigenvalues.

As the iteration goes, we will eventually converge to an upper triangular matrix form:

$\begin{split} A_k = R_kQ_k = \begin{bmatrix} \lambda_1 & X & \dots & X\\ 0 & \lambda_2 & \dots & X\\ & &\dots &\\ 0 & 0 & \dots & \lambda_n\\ \end{bmatrix}\end{split}$

where the diagonal values are the eigenvalues of the matrix. In each iteration of the QR method, factoring a matrix into an orthogonal and an upper triangular matrix can be done by using a special matrix called Householder matrix. We will not go into the mathematical details how you get the $$Q$$ and $$R$$ from the matrix, instead, we will use the Python function to obtain the two matrices directly.

TRY IT! Use the qr function in numpy.linalg to decompose matrix $$A = \begin{bmatrix} 0 & 2\\ 2 & 3\\ \end{bmatrix}$$. And verify the results.

import numpy as np
from numpy.linalg import qr

a = np.array([[0, 2],
[2, 3]])

q, r = qr(a)

print('Q:', q)
print('R:', r)

b = np.dot(q, r)
print('QR:', b)

Q: [[ 0. -1.]
[-1.  0.]]
R: [[-2. -3.]
[ 0. -2.]]
QR: [[0. 2.]
[2. 3.]]


TRY IT! Use the QR method to get the eigenvalues of matrix $$A = \begin{bmatrix} 0 & 2\\ 2 & 3\\ \end{bmatrix}$$. Do 20 iterations, and print out the 1st, 5th, 10th, and 20th iteration.

a = np.array([[0, 2],
[2, 3]])
p = [1, 5, 10, 20]
for i in range(20):
q, r = qr(a)
a = np.dot(r, q)
if i+1 in p:
print(f'Iteration {i+1}:')
print(a)

Iteration 1:
[[3. 2.]
[2. 0.]]
Iteration 5:
[[ 3.99998093  0.00976559]
[ 0.00976559 -0.99998093]]
Iteration 10:
[[ 4.00000000e+00  9.53674316e-06]
[ 9.53674316e-06 -1.00000000e+00]]
Iteration 20:
[[ 4.00000000e+00  9.09484250e-12]
[ 9.09494702e-12 -1.00000000e+00]]


We can see after the 5th iteration, the eigenvalues are converge to the correct ones.

In next section, we will take a look how to get the eigenvalues and eigenvectors in Python easily using the built-in function.