3.6.2 Diagonalisation and matrix powers

Working along the eigenvectors

Lesson 3.6.1 found the eigenvectors of

A = \begin{pmatrix} 1 & 2 \\ 3 & 2 \end{pmatrix}

They are $\mathbf{v}_1 = (2, 3)$ with $\lambda_1 = 4$ , and $\mathbf{v}_2 = (1, -1)$ with $\lambda_2 = -1$ . The two do not lie on one line, so they form a basis of the plane (lesson 3.5.1): every vector is a combination $c_1\mathbf{v}_1 + c_2\mathbf{v}_2$ . The numbers $c_1$ and $c_2$ are its coordinates in the eigenvector basis.

Take $\mathbf{x} = (4, 1)$ . The two equations $2c_1 + c_2 = 4$ and $3c_1 - c_2 = 1$ add up to $5c_1 = 5$ . So $c_1 = 1$ , $c_2 = 2$ , and $\mathbf{x} = \mathbf{v}_1 + 2\mathbf{v}_2$ .

A linear map keeps combinations (lesson 3.2.2), and it multiplies each eigenvector by its eigenvalue:

\begin{aligned} A\mathbf{x} &= A\mathbf{v}_1 + 2A\mathbf{v}_2 = 4\mathbf{v}_1 - 2\mathbf{v}_2 \\ &= (8, 12) + (-2, 2) = (6, 14) \end{aligned}

Multiplying directly agrees: $A(4, 1) = (4 + 2,\ 12 + 2) = (6, 14)$ . In eigenvector coordinates, $A$ turned $(1, 2)$ into $(4, -2)$ : it multiplied each coordinate by its own eigenvalue.

Each further multiplication does the same again, so $A^k\mathbf{x} = 4^k\mathbf{v}_1 + 2\,(-1)^k\,\mathbf{v}_2$ . Here $A^k$ means $A$ times itself $k$ times, and $A^0 = I$ .

$k$	0	1	2	3	4	5
$A^k\mathbf{x}$	$(4, 1)$	$(6, 14)$	$(34, 46)$	$(126, 194)$	$(514, 766)$	$(2046, 3074)$

For $k = 3$ the formula gives $64\,(2, 3) - 2\,(1, -1) = (126, 194)$ , without multiplying by $A$ three times.

Rule: Powers along the eigenvectors

\begin{gathered} \mathbf{x} = c_1\mathbf{v}_1 + c_2\mathbf{v}_2 \\[4pt] A^k\mathbf{x} = c_1\lambda_1^k\,\mathbf{v}_1 + c_2\lambda_2^k\,\mathbf{v}_2 \end{gathered}

Split $\mathbf{x}$ along the eigenvectors, raise each eigenvalue to the power $k$ , and add the parts back up. The same holds with $n$ terms for an $n \times n$ matrix with $n$ independent eigenvectors.

(c_1, c_2) = (1, 2)

import numpy as np A = np.array([[1, 2], [3, 2]]) v1, v2 = np.array([2, 3]), np.array([1, -1]) # eigenvalues 4 and -1 x = np.array([4, 1]) c = np.linalg.solve(np.column_stack([v1, v2]), x) print(c) for k in range(6): direct = np.linalg.matrix_power(A, k) @ x formula = c[0] * 4**k * v1 + c[1] * (-1)**k * v2 print(k, direct, formula)

⌘/Ctrl + Enter runs the cell

Diagonalisation

Put the eigenvectors side by side as the columns of a matrix $P$ . Put the eigenvalues, in the same order, on the diagonal of $\Lambda$ (capital lambda):

\begin{gathered} P = \begin{pmatrix} 2 & 1 \\ 3 & -1 \end{pmatrix} \\[6pt] \Lambda = \begin{pmatrix} 4 & 0 \\ 0 & -1 \end{pmatrix} \end{gathered}

Column $j$ of $AP$ is $A\mathbf{v}_j = \lambda_j\mathbf{v}_j$ , and so is column $j$ of $P\Lambda$ . So $AP = P\Lambda$ . When $P$ has an inverse, multiplying both sides by $P^{-1}$ on the right gives the factorisation.

Rule: Diagonalisation

A = P\Lambda P^{-1}

$P$ holds the eigenvectors as columns and $\Lambda$ the eigenvalues on its diagonal, in the same order. It works exactly when $A$ has $n$ independent eigenvectors, so that $P$ has an inverse; $A$ is then called diagonalisable. A matrix with $n$ different eigenvalues always is.

Read the product from right to left, as in lesson 3.2.2. $P^{-1}$ finds the coordinates of $\mathbf{x}$ along the eigenvectors, $\Lambda$ multiplies each coordinate by its eigenvalue, and $P$ builds the vector back.

$P$ has $ad - bc = -2 - 3 = -5$ , so the formula of lesson 3.2.3 gives

\begin{aligned} P^{-1} &= \frac{1}{-5}\begin{pmatrix} -1 & -1 \\ -3 & 2 \end{pmatrix} \\[4pt] &= \begin{pmatrix} 0.2 & 0.2 \\ 0.6 & -0.4 \end{pmatrix} \end{aligned}

and $P^{-1}(4, 1) = (1, 2)$ : the coordinates found above.

$P$ is not the only choice. Any nonzero multiples of the eigenvectors work, such as $(1, 1.5)$ for $(2, 3)$ , and so does another order, as long as $\Lambda$ follows it.

AP = P\Lambda

import numpy as np A = np.array([[1, 2], [3, 2]]) P = np.array([[2, 1], [3, -1]]) # eigenvectors (2, 3) and (1, -1) as columns Lam = np.diag([4, -1]) # eigenvalues in the same order P_inv = np.linalg.inv(P) print(np.round(P_inv, 4)) print(np.array_equal(A @ P, P @ Lam), np.allclose(P @ Lam @ P_inv, A)) print(np.round(P_inv @ np.array([4, 1]), 6))

⌘/Ctrl + Enter runs the cell

Matrix powers

Squaring the factorisation, the $P^{-1}P$ in the middle is $I$ and drops out:

\begin{aligned} A^2 &= P\Lambda P^{-1}\,P\Lambda P^{-1} \\ &= P\Lambda^2 P^{-1} \end{aligned}

Every further factor does the same, so only $\Lambda$ is raised to the power. A power of a diagonal matrix is just the powers of its diagonal entries.

Rule: Powers of a diagonalisable matrix

\begin{gathered} A^k = P\Lambda^k P^{-1} \\[4pt] \Lambda^k = \begin{pmatrix} \lambda_1^k & 0 \\ 0 & \lambda_2^k \end{pmatrix} \end{gathered}

$A^k$ has the same eigenvectors as $A$ , and eigenvalues $\lambda_1^k$ and $\lambda_2^k$ . If no eigenvalue is 0, the inverse works the same way: $A^{-1} = P\Lambda^{-1}P^{-1}$ .

For $k = 5$ , $\Lambda^5$ has the diagonal $4^5 = 1024$ and $(-1)^5 = -1$ :

\begin{gathered} P\Lambda^5 = \begin{pmatrix} 2048 & -1 \\ 3072 & 1 \end{pmatrix} \\[6pt] A^5 = P\Lambda^5 P^{-1} = \begin{pmatrix} 409 & 410 \\ 615 & 614 \end{pmatrix} \end{gathered}

For $k = 10$ the work is the same: two numbers to raise, two products to form. The answer is $\begin{pmatrix} 419431 & 419430 \\ 629145 & 629146 \end{pmatrix}$ , which nine multiplications by $A$ confirm.

The inverse has eigenvalues $\frac{1}{4}$ and $-1$ , and $P\Lambda^{-1}P^{-1} = \begin{pmatrix} -0.5 & 0.5 \\ 0.75 & -0.25 \end{pmatrix}$ , the same as np.linalg.inv(A).

k = 5

import numpy as np A = np.array([[1, 2], [3, 2]]) P = np.array([[2, 1], [3, -1]]) P_inv = np.linalg.inv(P) for k in [5, 10]: Lam_k = np.diag([4.0**k, (-1.0)**k]) print(np.round(P @ Lam_k @ P_inv, 6) + 0.0) print(np.linalg.matrix_power(A, k)) print(np.round(P @ np.diag([1 / 4, -1.0]) @ P_inv, 6), np.round(np.linalg.inv(A), 6))

⌘/Ctrl + Enter runs the cell

The largest eigenvalue takes over

In $A^k\mathbf{x} = 4^k\mathbf{v}_1 + 2\,(-1)^k\,\mathbf{v}_2$ , the first part grows four times at each step while the second keeps its size. Divide by $4^k$ :

$k$	0	1	2	3	4	5
$A^k\mathbf{x} / 4^k$	$(4, 1)$	$(1.5, 3.5)$	$(2.125, 2.875)$	$(1.9688, 3.0312)$	$(2.0078, 2.9922)$	$(1.998, 3.002)$

The rows close in on $\mathbf{v}_1 = (2, 3)$ , from alternate sides because $\lambda_2$ is negative. What is left of $\mathbf{v}_2$ is $2\left(-\frac{1}{4}\right)^k\mathbf{v}_2$ , a quarter of its size at each step.

Rule: The dominant eigenvalue

\begin{gathered} A^k\mathbf{x} \approx c_1\lambda_1^k\,\mathbf{v}_1 \\[4pt] \text{when } \lvert\lambda_1\rvert > \lvert\lambda_2\rvert \end{gathered}

For large $k$ , the eigenvalue of largest size wins, as long as $c_1 \neq 0$ . The rest shrinks like $\lvert \lambda_2 / \lambda_1 \rvert^k$ . In the long run $A^k\mathbf{x}$ grows if $\lvert\lambda_1\rvert > 1$ , settles if $\lambda_1 = 1$ , and dies away if $\lvert\lambda_1\rvert < 1$ .

A^k\mathbf{x}

import numpy as np A = np.array([[1, 2], [3, 2]]) x = np.array([4, 1]) v1, v2 = np.array([2, 3]), np.array([1, -1]) for k in range(6): scaled = np.linalg.matrix_power(A, k) @ x / 4**k print(k, np.round(scaled, 4), np.round(scaled - v1, 4) + 0.0)

⌘/Ctrl + Enter runs the cell

Two towns that swap people

Each year 10% of the people in town A move to town B, and 20% of those in town B move to town A. Write $\mathbf{s} = (s_A, s_B)$ for the shares of the people living in each town, so $s_A + s_B = 1$ .

Next year A keeps 90% of its own and gains 20% of B's: $0.9s_A + 0.2s_B$ . B gets the rest: $0.1s_A + 0.8s_B$ . So each year multiplies $\mathbf{s}$ by

M = \begin{pmatrix} 0.9 & 0.2 \\ 0.1 & 0.8 \end{pmatrix}

Column $j$ says where the people of town $j$ go. Each column adds up to 1, so nobody is lost and the shares still add up to 1.

$M$ has trace 1.7 and determinant $0.72 - 0.02 = 0.7$ . Its polynomial is $\lambda^2 - 1.7\lambda + 0.7 = (\lambda - 1)(\lambda - 0.7)$ , so the eigenvalues are 1 and 0.7.

$\lambda_1 = 1$ : the shortcut of lesson 3.6.1 gives $(0.2, 0.1)$ , or as shares, $\left(\frac{2}{3}, \frac{1}{3}\right)$ . Here $M\mathbf{s} = \mathbf{s}$ : a steady state that stays put year after year.
$\lambda_2 = 0.7$ : the eigenvector $(1, -1)$ moves people from one town to the other. It shrinks by 0.7 each year.

Start from $\mathbf{s}_0 = (0.5, 0.5)$ , which is $\left(\frac{2}{3}, \frac{1}{3}\right) - \frac{1}{6}(1, -1)$ . After $k$ years the shares are

\mathbf{s}_k = \left(\tfrac{2}{3}, \tfrac{1}{3}\right) - \tfrac{1}{6}\,(0.7)^k\,(1, -1)

After 10 years, $0.7^{10} \approx 0.0282$ and $\mathbf{s}_{10} \approx (0.662, 0.338)$ . The gap to the steady state is multiplied by 0.7 each year, so it halves about every two years, since $0.7^2 = 0.49$ .

From any start, the share in town A settles on $\frac{2}{3}$ , the eigenvector of $\lambda = 1$ . On a log scale the gaps are straight lines: each year multiplies them by the other eigenvalue, 0.7.

Lesson 5.11.2 reads such shares as probabilities: $M$ becomes a Markov chain, and the eigenvector of $\lambda = 1$ its long-run distribution.

\lambda = 1

import numpy as np M = np.array([[0.9, 0.2], [0.1, 0.8]]) # column j: where town j's people go s = np.array([0.5, 0.5]) for year in range(11): if year in (0, 1, 2, 3, 5, 10): print(year, np.round(s, 4)) s = M @ s values, V = np.linalg.eig(M) values, V = values.real, V.real order = np.argsort(values)[::-1] values, V = values[order], V[:, order] steady = V[:, 0] / V[:, 0].sum() # the lambda = 1 eigenvector, as shares print(np.round(values, 6), np.round(steady, 4)) gap = 0.5 - 2 / 3 print(np.round(np.array([2 / 3, 1 / 3]) + gap * 0.7**10 * np.array([1, -1]), 4))

⌘/Ctrl + Enter runs the cell

When a matrix cannot be diagonalised

The shear $H = \begin{pmatrix} 1 & 1 \\ 0 & 1 \end{pmatrix}$ has the eigenvalue 1 twice, and only one eigenvector line, the $x_1$ -axis (lesson 3.6.1). Any two of its eigenvectors lie on that line, so a matrix $P$ built from them is singular. $H$ is not diagonalisable.

Its powers show it. Three shears make a shear by 3 (lesson 3.2.2), and in general $H^k = \begin{pmatrix} 1 & k \\ 0 & 1 \end{pmatrix}$ . The corner grows by 1 at every step, which no power of the eigenvalue 1 can produce.

A repeated eigenvalue alone is not the problem. $2I$ also has one eigenvalue twice, and it is already diagonal. What matters is whether there are $n$ independent eigenvectors.

Complex eigenvalues are another matter. The quarter turn $R$ has eigenvalues $i$ and $-i$ , and it can be diagonalised with complex numbers, which is what eig does. For example $R^4 = I$ , matching $i^4 = 1$ (lesson 1.7.2).

H

import numpy as np H = np.array([[1, 1], [0, 1]]) # the shear values, V = np.linalg.eig(H) V = V.real * np.sign(V.real[0]) print(np.round(values.real, 6), np.round(V, 4) + 0.0) print(abs(np.linalg.det(V)) < 1e-12) # P would be singular for k in [2, 3, 10]: print(k, np.linalg.matrix_power(H, k).tolist()) R = np.array([[0, -1], [1, 0]]) # the quarter turn print(np.linalg.matrix_power(R, 4))

⌘/Ctrl + Enter runs the cell

The power method

The dominant eigenvalue gives a way to find it without the characteristic polynomial. Multiply a starting vector by $A$ again and again: its direction settles on $\mathbf{v}_1$ . Rescale it to length 1 at every step, so the numbers neither explode nor vanish.

To read off $\lambda_1$ , use the dot product $\mathbf{x} \cdot (A\mathbf{x})$ for the current unit vector $\mathbf{x}$ . If $\mathbf{x}$ were exactly a unit eigenvector, this would be $\mathbf{x} \cdot \lambda_1\mathbf{x} = \lambda_1$ .

Rule: The power method

\begin{gathered} \mathbf{x} \leftarrow \frac{A\mathbf{x}}{\lVert A\mathbf{x} \rVert} \\[4pt] \lambda_1 \approx \mathbf{x} \cdot (A\mathbf{x}) \end{gathered}

Repeat until the estimate stops changing. The error shrinks like $\lvert \lambda_2 / \lambda_1 \rvert^k$ , so it is fast when $\lambda_1$ is far ahead and slow when $\lambda_2$ is close behind.

From $(1, 0)$ , the estimates for $A$ run 3.4, 4.0462, 3.9804, 4.0044, 3.9989, … The error falls about four times at each step, since $\lvert \lambda_2/\lambda_1 \rvert = \frac{1}{4}$ . After 8 steps $\mathbf{x} \approx (0.5547, 0.832)$ , which is $(2, 3)$ scaled to length 1.

The method fails in two cases:

No single winner. If $\lvert\lambda_1\rvert = \lvert\lambda_2\rvert$ , as for the reflection $\begin{pmatrix} 1 & 0 \\ 0 & -1 \end{pmatrix}$ or any turn, the direction never settles.
A start with $c_1 = 0$ . From $(1, -1)$ , $A$ only flips the vector, and the method finds $\lambda_2 = -1$ . A random start almost never has $c_1 = 0$ exactly.

NumPy's eig uses a refined descendant of this idea, the QR algorithm. The power method itself is used for very large matrices when only the top eigenvalue is needed. Practice set 3.6.4 asks you to write it.

A

import numpy as np A = np.array([[1, 2], [3, 2]]) x = np.array([1.0, 0.0]) for step in range(1, 9): y = A @ x x = y / np.linalg.norm(y) # rescale to length 1 estimate = x @ (A @ x) print(step, np.round(x, 4), round(estimate, 4), f"{abs(estimate - 4):.1e}")

⌘/Ctrl + Enter runs the cell

Common mistakes

Match the order. Column $j$ of $P$ must be the eigenvector of the $j$ -th diagonal entry of $\Lambda$ .
$P\Lambda P^{-1}$ , not $P^{-1}\Lambda P$ . Reading right to left, $P^{-1}$ acts first and finds the eigenvector coordinates.
A matrix power is not a power of each entry. A**3 in NumPy cubes every entry. Use np.linalg.matrix_power(A, 3) or A @ A @ A.
Size decides who wins. With eigenvalues $-5$ and 3, the $-5$ takes over, and $A^k\mathbf{x}$ flips direction at every step.
A repeated eigenvalue does not settle the question. $2I$ is diagonalisable and the shear is not: count the independent eigenvectors.
The power method needs a clear winner and a start with some part along $\mathbf{v}_1$ .

Diagonalisation and matrix powers

0 of 12 correct

1.
$A\mathbf{v} = 3\mathbf{v}$ for $\mathbf{v} = (1, 2)$ . Find the second entry of $A^4\mathbf{v}$
2.
$\mathbf{x} = 2\mathbf{v}_1 + \mathbf{v}_2$ , with $A\mathbf{v}_1 = 3\mathbf{v}_1$ and $A\mathbf{v}_2 = -\mathbf{v}_2$ . Then $A^2\mathbf{x} = a\mathbf{v}_1 + b\mathbf{v}_2$ . Find $a$
3.
The same $\mathbf{x}$ : find $b$ in $A^3\mathbf{x} = a\mathbf{v}_1 + b\mathbf{v}_2$
4.
$A = P\Lambda P^{-1}$ with $\Lambda = \begin{pmatrix} 2 & 0 \\ 0 & 0.5 \end{pmatrix}$ . What is the top-left entry of $\Lambda^6$ ?
5.
$A = \begin{pmatrix} 1 & 2 \\ 3 & 2 \end{pmatrix}$ has eigenvalues 4 and $-1$ . Find the trace of $A^5$
6.
A $2 \times 2$ matrix has eigenvalues 2 and 5. Find $\det(A^3)$
7.
Each year 20% of town A's people move to B and 30% of B's move to A, so $M = \begin{pmatrix} 0.8 & 0.3 \\ 0.2 & 0.7 \end{pmatrix}$ . What is the steady share of town A?
8.
For the same $M$ , by what factor does the gap to the steady state shrink each year?
9.
$\mathbf{x} = (0.6, 0.8)$ has length 1, and $A\mathbf{x} = (2.4, 3.2)$ . Use the power method's estimate to find the eigenvalue
10.
The eigenvalues of $A$ are 5 and 1. By what factor does the power method's error shrink at each step?
11.
Find the top-right entry of $\begin{pmatrix} 1 & 1 \\ 0 & 1 \end{pmatrix}^{20}$
12.
A $3 \times 3$ matrix has eigenvalues 1, 2 and 3. Is it diagonalisable? Type 1 for yes, 0 for no

How to type answers

Powers with ^: x^7, 2^(1/3), a^-2.
Multiply with * or by writing side by side: 3x^2, 2(x+1).
ln(x) is the natural log, log(x, 2) is log base 2, e is Euler's number, sqrt(x) is a square root.
sin(x), cos(x) and tan(x) take radians; write pi for π and sin(x)^2 for sin²x.
Your first answer to each question is saved: it counts towards your score and brings the question back for review.
Exact answers like ln(20)/ln(3) are best. Decimals need about 4 significant figures.

\mathbf{x} = c_1\mathbf{v}_1 + c_2\mathbf{v}_2

Diagonalisation and matrix powers

Working along the eigenvectors

Rule: Powers along the eigenvectors

Diagonalisation

Rule: Diagonalisation

Matrix powers

Rule: Powers of a diagonalisable matrix

The largest eigenvalue takes over

Rule: The dominant eigenvalue

Two towns that swap people

When a matrix cannot be diagonalised

The power method

Rule: The power method

Common mistakes

Diagonalisation and matrix powers

Summary card

Free further reading