# Chapter 94. QR Algorithm

# Chapter 94. QR Algorithm

The QR algorithm is an iterative method for computing eigenvalues of a square matrix. It is one of the central algorithms of numerical linear algebra.

The method is based on a simple idea. Factor the current matrix into an orthogonal factor and an upper triangular factor,

$$
A_k = Q_kR_k,
$$

then reverse the factors:

$$
A_{k+1}=R_kQ_k.
$$

Since

$$
A_{k+1} =
R_kQ_k =
Q_k^TA_kQ_k,
$$

the new matrix is orthogonally similar to the old one. Therefore it has the same eigenvalues. Repeating this process transforms the matrix toward a form where the eigenvalues can be read from the diagonal or from small diagonal blocks. The practical QR algorithm is usually implemented after reduction to Hessenberg form and with shifts and deflation.

## 94.1 The Eigenvalue Problem

Let

$$
A\in \mathbb{R}^{n\times n}.
$$

The eigenvalue problem asks for scalars \(\lambda\) and nonzero vectors \(v\) such that

$$
Av=\lambda v.
$$

Equivalently, the eigenvalues are the roots of

$$
\det(A-\lambda I)=0.
$$

For small matrices, this equation may be solved directly. For large matrices, computing and solving the characteristic polynomial is numerically poor. The QR algorithm avoids explicit polynomial root-finding.

Its goal is to transform \(A\) by similarity transformations into a nearly triangular matrix. Since triangular matrices have eigenvalues on the diagonal, the eigenvalues then become visible.

## 94.2 Similarity Transformations

Two matrices \(A\) and \(B\) are similar if

$$
B=S^{-1}AS
$$

for some invertible matrix \(S\).

Similar matrices have the same eigenvalues.

The QR algorithm uses orthogonal similarity transformations:

$$
A_{k+1}=Q_k^TA_kQ_k.
$$

This is numerically important because orthogonal transformations preserve vector norms and are stable in floating point arithmetic.

The algorithm therefore changes the representation of the matrix without changing its eigenvalues.

## 94.3 QR Factorization Step

At iteration \(k\), compute a QR factorization

$$
A_k=Q_kR_k,
$$

where:

| Symbol | Meaning |
|---|---|
| \(Q_k\) | orthogonal matrix |
| \(R_k\) | upper triangular matrix |

Then define

$$
A_{k+1}=R_kQ_k.
$$

Because

$$
R_k=Q_k^TA_k,
$$

we get

$$
A_{k+1}=Q_k^TA_kQ_k.
$$

Thus every QR step is an orthogonal similarity transformation.

The sequence

$$
A_0,A_1,A_2,\ldots
$$

preserves all eigenvalues of the original matrix.

## 94.4 Basic QR Algorithm

The unshifted QR algorithm is:

```text id="v21zro"
qr_algorithm(A, max_iter):
    A_k = A

    for k = 0 to max_iter - 1:
        Q, R = qr_factorization(A_k)
        A_k = R * Q

    return A_k
```

If the method converges, \(A_k\) approaches an upper triangular matrix:

$$
A_k \to T.
$$

The eigenvalues are then approximately

$$
t_{11},t_{22},\ldots,t_{nn}.
$$

For symmetric matrices, the limit is diagonal rather than merely triangular.

## 94.5 Why Reversing the Factors Works

The QR step starts from

$$
A_k=Q_kR_k.
$$

Multiplying on the left by \(Q_k^T\) gives

$$
R_k=Q_k^TA_k.
$$

Then

$$
A_{k+1}=R_kQ_k =
Q_k^TA_kQ_k.
$$

Thus the operation \(A_k\mapsto A_{k+1}\) is not an arbitrary multiplication. It is a change of orthonormal basis.

The eigenvalues remain fixed, while the matrix representation changes.

The iteration attempts to choose these basis changes so that the matrix becomes closer to triangular form.

## 94.6 Schur Form

The Schur decomposition states that a square complex matrix can be written as

$$
A=QTQ^*,
$$

where \(Q\) is unitary and \(T\) is upper triangular.

For real matrices, the real Schur form is quasi-upper triangular. It has \(1\times 1\) blocks for real eigenvalues and \(2\times 2\) blocks for complex conjugate pairs.

The QR algorithm may be understood as a method for computing Schur form. Once Schur form is obtained, eigenvalues are read from the diagonal entries or from the \(2\times 2\) diagonal blocks.

## 94.7 Symmetric Case

If

$$
A=A^T,
$$

then all eigenvalues are real, and the Schur form is diagonal.

The QR algorithm attempts to compute

$$
A=Q\Lambda Q^T,
$$

where

$$
\Lambda=
\operatorname{diag}(\lambda_1,\lambda_2,\ldots,\lambda_n).
$$

In this case, the iterates tend toward a diagonal matrix under suitable conditions.

The diagonal entries approximate the eigenvalues. The accumulated orthogonal transformations approximate the eigenvectors.

## 94.8 Relation to Power Iteration

Power iteration repeatedly applies \(A\) to one vector and normalizes:

$$
x,\quad Ax,\quad A^2x,\ldots.
$$

The QR algorithm may be viewed as a simultaneous version of power iteration.

Instead of tracking one vector, it tracks a full orthonormal basis.

Repeated multiplication by \(A\) tends to align basis vectors with invariant directions. QR factorization reorthogonalizes the basis at each step.

This explains the connection between eigenvalue separation and convergence. When eigenvalues are well separated, directions separate faster.

## 94.9 Accumulating Eigenvectors

If eigenvectors are desired, one accumulates the orthogonal transformations.

Let

$$
Z_0=I.
$$

At each step,

$$
A_k=Q_kR_k,
\qquad
A_{k+1}=R_kQ_k,
$$

and update

$$
Z_{k+1}=Z_kQ_k.
$$

Then

$$
A_k=Z_k^TAZ_k.
$$

If \(A_k\) converges to a diagonal matrix \(\Lambda\), then

$$
A\approx Z_k\Lambda Z_k^T.
$$

The columns of \(Z_k\) approximate eigenvectors in the symmetric case.

## 94.10 Cost of the Naive Algorithm

A direct QR factorization of a dense \(n\times n\) matrix costs

$$
O(n^3)
$$

operations per iteration.

Many iterations may be required.

Therefore the naive algorithm is too expensive for practical dense eigenvalue computation.

The practical QR algorithm first reduces the matrix to a simpler form:

| Matrix type | Reduced form |
|---|---|
| General square matrix | upper Hessenberg form |
| Symmetric matrix | tridiagonal form |

The QR iteration then preserves this reduced structure and becomes much cheaper.

## 94.11 Hessenberg Form

An upper Hessenberg matrix has zeros below the first subdiagonal:

$$
H=
\begin{bmatrix}
h_{11} & h_{12} & h_{13} & \cdots & h_{1n}\\
h_{21} & h_{22} & h_{23} & \cdots & h_{2n}\\
0 & h_{32} & h_{33} & \cdots & h_{3n}\\
0 & 0 & h_{43} & \cdots & h_{4n}\\
\vdots & \vdots & \vdots & \ddots & \vdots\\
0 & 0 & 0 & h_{n,n-1} & h_{nn}
\end{bmatrix}.
$$

Every square matrix can be reduced to Hessenberg form by orthogonal similarity transformations.

For a general matrix, practical QR iteration is applied to this Hessenberg matrix rather than to the original dense matrix. This reduction is a standard first stage of the QR algorithm.

## 94.12 Tridiagonal Form for Symmetric Matrices

If \(A\) is symmetric, Hessenberg reduction preserves symmetry.

A symmetric Hessenberg matrix is tridiagonal:

$$
T=
\begin{bmatrix}
d_1 & e_1 & 0 & \cdots & 0\\
e_1 & d_2 & e_2 & \cdots & 0\\
0 & e_2 & d_3 & \cdots & 0\\
\vdots & \vdots & \vdots & \ddots & e_{n-1}\\
0 & 0 & 0 & e_{n-1} & d_n
\end{bmatrix}.
$$

Thus the symmetric QR algorithm usually has two stages:

| Stage | Operation |
|---|---|
| Reduction | \(A\) is reduced to tridiagonal form |
| Iteration | QR steps are applied to the tridiagonal matrix |

This structure makes the symmetric eigenvalue problem much cheaper than the fully general case.

## 94.13 Preservation of Hessenberg Structure

If \(H_k\) is upper Hessenberg, then a QR step produces another upper Hessenberg matrix:

$$
H_{k+1}=R_kQ_k.
$$

This is crucial.

It means the expensive reduction is performed once. Later iterations stay within the Hessenberg class.

For a Hessenberg matrix, one QR step can be computed in

$$
O(n^2)
$$

operations.

For a symmetric tridiagonal matrix, one step can be computed in

$$
O(n)
$$

operations.

## 94.14 Shifted QR Algorithm

Unshifted QR iteration may converge slowly.

A shifted QR step chooses a scalar \(\mu_k\), factors

$$
A_k-\mu_k I=Q_kR_k,
$$

then defines

$$
A_{k+1}=R_kQ_k+\mu_k I.
$$

This is still a similarity transformation:

$$
A_{k+1} =
Q_k^TA_kQ_k.
$$

The shift changes the factorization step but preserves eigenvalues.

The purpose of \(\mu_k\) is to accelerate convergence by choosing a value close to an eigenvalue. Practical QR algorithms use shifts extensively.

## 94.15 Shifted QR Pseudocode

```text id="dbv8yz"
shifted_qr_algorithm(A, max_iter):
    A_k = hessenberg_reduction(A)

    for k = 0 to max_iter - 1:
        mu = choose_shift(A_k)

        Q, R = qr_factorization(A_k - mu * I)
        A_k = R * Q + mu * I

        deflate if a subdiagonal entry is small

    return A_k
```

This pseudocode hides important implementation details. Modern algorithms usually use implicit shifts rather than explicitly forming a full QR factorization at every step.

## 94.16 Rayleigh Quotient Shift

A simple shift choice is

$$
\mu_k=(A_k)_{nn}.
$$

This uses the bottom-right diagonal entry as an eigenvalue estimate.

For some matrices, this accelerates convergence substantially.

However, more refined shifts are usually preferred.

## 94.17 Wilkinson Shift

For symmetric tridiagonal matrices, the Wilkinson shift is a standard choice.

It uses the eigenvalue of the trailing \(2\times 2\) principal submatrix that is closer to the bottom-right diagonal entry.

If the trailing block is

$$
\begin{bmatrix}
a & b\\
b & d
\end{bmatrix},
$$

then the shift is chosen from its two eigenvalues.

This often gives very rapid convergence of the bottom subdiagonal entry to zero.

## 94.18 Deflation

Deflation occurs when a subdiagonal entry becomes sufficiently small.

Suppose

$$
|(A_k)_{i+1,i}|
$$

is tiny.

Then \(A_k\) is nearly block upper triangular:

$$
A_k
\approx
\begin{bmatrix}
A_{11} & A_{12}\\
0 & A_{22}
\end{bmatrix}.
$$

The eigenvalues of \(A_k\) are then approximately the union of the eigenvalues of \(A_{11}\) and \(A_{22}\).

The problem can be split into smaller eigenvalue problems.

In practice, deflation is essential for efficiency. It allows one eigenvalue or one conjugate pair to be isolated and removed from the active iteration.

## 94.19 Practical Deflation Test

A common deflation test checks whether

$$
|a_{i+1,i}|
\le
\tau \left(|a_{ii}|+|a_{i+1,i+1}|\right),
$$

where \(\tau\) is a tolerance related to machine precision.

If this inequality holds, the subdiagonal entry is set to zero.

This produces a block split.

The QR algorithm is then applied separately to each active block.

## 94.20 Real Schur Form and Complex Eigenvalues

A real matrix may have complex eigenvalues.

Since the QR algorithm can be performed using real arithmetic, complex conjugate eigenvalues appear as \(2\times 2\) diagonal blocks in real Schur form.

A block

$$
\begin{bmatrix}
a & b\\
c & d
\end{bmatrix}
$$

represents the eigenvalues of that \(2\times 2\) matrix.

If its characteristic polynomial has negative discriminant, the corresponding eigenvalues are complex conjugates.

Thus real QR iteration does not need to store complex numbers to represent complex eigenvalue pairs.

## 94.21 Implicit QR Algorithm

Modern implementations usually use the implicit QR algorithm.

Instead of explicitly computing

$$
Q_kR_k=A_k-\mu I
$$

and multiplying

$$
R_kQ_k+\mu I,
$$

the algorithm applies a sequence of orthogonal transformations directly to the Hessenberg matrix.

These transformations introduce a small nonzero pattern called a bulge and then chase it down the matrix until Hessenberg form is restored.

This process is called bulge chasing.

Implicit QR is more efficient, more stable, and better suited to multiple shifts.

## 94.22 Francis QR Step

The implicit shifted QR step is often called the Francis QR step.

It was developed independently by John Francis and Vera Kublanovskaya in the late 1950s and early 1960s. The QR algorithm is historically one of the major developments in numerical eigenvalue computation.

In modern dense eigenvalue solvers, the Francis step is the practical core of the QR algorithm.

## 94.23 Double Shifts

For real matrices, complex eigenvalues occur in conjugate pairs.

A double-shift QR step uses a quadratic polynomial

$$
p(t)=(t-\mu)(t-\overline{\mu})
$$

so that the computation can remain real even when the shifts are complex.

The shifts are often chosen from the trailing \(2\times 2\) block.

Double shifts are standard in real nonsymmetric QR algorithms.

## 94.24 Numerical Stability

The QR algorithm uses orthogonal transformations.

Orthogonal transformations preserve norms:

$$
\|Qx\|_2=\|x\|_2.
$$

This gives the method strong numerical stability compared with methods based on nonorthogonal transformations.

The LR algorithm, an older method based on LU factorization, is generally less stable. QR replaced it in practical eigenvalue computation largely because orthogonal transformations behave better numerically.

## 94.25 Eigenvectors

The QR algorithm naturally computes eigenvalues.

Eigenvectors require additional work.

For symmetric matrices, if the accumulated orthogonal matrix \(Z\) is stored, then

$$
A\approx Z\Lambda Z^T,
$$

and the columns of \(Z\) approximate eigenvectors.

For nonsymmetric matrices, eigenvectors can be computed from Schur form by solving triangular systems.

Often numerical software first computes Schur form, then optionally computes eigenvectors.

## 94.26 Relation to Hessenberg and Schur Decompositions

The practical dense eigenvalue pipeline is:

| Step | Purpose |
|---|---|
| Balance, optional | improve scaling |
| Hessenberg reduction | reduce cost of iteration |
| Implicit shifted QR | converge toward Schur form |
| Deflation | split completed blocks |
| Back transformation | recover Schur vectors or eigenvectors |

For symmetric matrices, the corresponding pipeline is:

| Step | Purpose |
|---|---|
| Tridiagonal reduction | reduce dense symmetric problem |
| Shifted QR or related method | compute eigenvalues |
| Accumulate transformations | compute eigenvectors |

The QR algorithm is therefore part of a larger eigenvalue solver.

## 94.27 Comparison with Power Iteration

| Feature | Power iteration | QR algorithm |
|---|---|---|
| Goal | dominant eigenpair | full spectrum |
| Main operation | matrix-vector product | orthogonal similarity steps |
| Memory | low | higher |
| Best use | large sparse leading eigenvalue | dense full eigenvalue problem |
| Convergence | depends on dominant gap | improved by shifts and deflation |
| Eigenvectors | one or few | many or all |

Power iteration is simple and sparse-friendly. QR is more complete and is the standard dense eigenvalue method.

## 94.28 Comparison with Krylov Methods

Krylov methods such as Arnoldi and Lanczos are preferred for large sparse matrices when only a few eigenvalues are needed.

The QR algorithm is preferred for dense matrices when most or all eigenvalues are needed.

| Setting | Typical choice |
|---|---|
| Dense matrix, all eigenvalues | QR algorithm |
| Symmetric dense matrix | tridiagonal reduction plus QR or related methods |
| Sparse matrix, few eigenvalues | Lanczos or Arnoldi |
| Matrix-free operator | Krylov method |
| Full Schur form | QR algorithm |

Thus QR and Krylov methods serve different computational regimes.

## 94.29 Failure Modes and Limitations

The QR algorithm is robust, but practical implementations must handle several issues.

| Issue | Handling |
|---|---|
| Slow convergence | use shifts |
| Complex eigenvalues | use real Schur \(2\times 2\) blocks |
| High per-step cost | reduce to Hessenberg or tridiagonal form |
| Completed eigenvalues | deflate |
| Roundoff | use orthogonal transformations |
| Poor scaling | balance matrix before iteration |

The naive unshifted QR algorithm is mainly pedagogical. Practical QR algorithms use Hessenberg reduction, implicit shifts, and deflation.

## 94.30 Summary

The QR algorithm computes eigenvalues by repeatedly applying QR factorization and reversing the factors:

$$
A_k=Q_kR_k,
\qquad
A_{k+1}=R_kQ_k.
$$

Each step is an orthogonal similarity transformation:

$$
A_{k+1}=Q_k^TA_kQ_k.
$$

Therefore the eigenvalues are preserved.

Under suitable conditions, the iterates approach Schur form. For symmetric matrices, they approach diagonal form. In practical implementations, the algorithm is made efficient by:

| Technique | Role |
|---|---|
| Hessenberg reduction | lowers cost for general matrices |
| Tridiagonal reduction | lowers cost for symmetric matrices |
| Shifts | accelerate convergence |
| Deflation | splits completed blocks |
| Implicit QR | avoids explicit expensive steps |
| Orthogonal transformations | provide numerical stability |

The QR algorithm is the classical dense-matrix method for computing the full eigenvalue spectrum.