Householder Tridiagonalization

M.P.
2022-05-11

Introduction

Useful matrix manipulations through Householder Reflections

(Photo Source: Townsend)

Invented by American Mathematician A.S. Householder in 1958 (Kreyszig).

Orthogonal projections generate Householder matrices, \( P \), which zero-out values of \( A \) when multiplied (Orthogonal Matrices: Lecture 8).

Introduction to Householder Tridiagonalization

Technique to tridiagonalize a symmetric matrix (Kreyszig).
Tridiagonal matrix : nonzero entries on main and adjacent diagonals, efficient for storage (Weisstein).

\[ T=\begin{bmatrix} 5&5&0&0\\4&1&3&0\\0&9&7&10\\0&0&-1&-5 \end{bmatrix} \]

\( T \) is similar to \( A \): eigenvalues unchanged and there exist \( S,S^{-1} \) such that \( T = SAS^{-1} \) (Rabinoff & Margalit).

Description of Method

From an \( n\times n \) symmetric matrix \( A \), must create \( n-2 \) Householder matrices, \( P \):

\[ \begin{equation}P = I - 2\mathbf{u}\mathbf{u}^T\end{equation} \]

\( P \) are symmetric (\( P^T = P \)) and orthogonal ( \( P^T=P^{-1} \))(Lee)
Compute matrix product in each of \( j \) steps to update \( A \):

\[ A = P_jAP_j^T \]

After \( n-2 \) multiplications, the matrix \( A \) will be tridiagonal (Kreyszig).

Description of Method continued ...

\( P_j = I-2\mathbf{u}_j\mathbf{u}_j^T, 1\leq j\leq n-2 \)

\( \mathbf{u}_j \) constructed from corresponding columns of A, \( \mathbf{a}_j \).
First \( j \) entries of \( \mathbf{u}_j \) are zero, rest are calculated from formulas (see Kreyszig 20.9):

\[ \begin{array}{l} &u_{j+1} = \sqrt{\frac{1}{2}\left( 1+\frac{|a_{j+1,j}|}{s_j}\right)}\\ &u_i = \displaystyle\frac{a_{i,j} \ \mathrm{sgn} (a_{j+1,j})}{2u_{j+1}s_j}, \ j+2 \leq i\leq n \end{array} \]

\[ s_j = \sqrt{ \sum_{i = j+1}^n(a_{i,j})^2} \]

\[ \mathbf{u} = \begin{bmatrix} 0\\ 0\\ \vdots \\ 0\\ u_{j+1}\\ \vdots\\ u_n \end{bmatrix} \]

Octave Code

Declare symmetric matrix A
A =

   4   3   2   1
   3   2   1   4
   2   1   4   3
   1   4   3   2

Vector u calculated from each column of A with specific formulas
u =

        0
   0.9492
   0.2816
   0.1408

Octave Code continued...

Householder matrix P = I-2uu': 
P =

   1.0000        0        0        0
        0  -0.8018  -0.5345  -0.2673
        0  -0.5345   0.8414  -0.0793
        0  -0.2673  -0.0793   0.9604

Left multiply A by P and right-multiply by P': 
A =

   4.0000  -3.7417  -0.0000        0
  -3.7417   6.0000  -1.0102  -3.9795
  -0.0000  -1.0102   2.4552   0.1586
        0  -3.9795   0.1586  -0.4552

Octave Code continued...

A new u vector is calculated from the next column of A: 
u =

        0
        0
   0.7893
   0.6140

New Householder matrix P = I-2uu': 
P =

   1.0000        0        0        0
        0   1.0000        0        0
        0        0  -0.2461  -0.9693
        0        0  -0.9693   0.2461

Octave Code continued ...

Updated matrix A is calculated from PAP'
A =

   4.0000e+00  -3.7417e+00   4.7806e-17   1.8832e-16
  -3.7417e+00   6.0000e+00   4.1057e+00   1.3323e-15
   4.7806e-17   4.1057e+00  -2.0339e-01   8.3349e-01
   1.8832e-16   9.9920e-16   8.3349e-01   2.2034e+00

Commands and Output

Values off main diagonals on the order of 10^-15, so we round to zero:
T =

   4.0000  -3.7417        0        0
  -3.7417   6.0000   4.1057        0
        0   4.1057  -0.2034   0.8335
        0        0   0.8335   2.2034

Original matrix
A =

   4   3   2   1
   3   2   1   4
   2   1   4   3
   1   4   3   2

Commands and Output continued...

To prove that we get a similar matrix, use Octave's eig() function to calculate the eigenvalues and compare.
eigT =

   10.0000
   -2.8284
    2.0000
    2.8284

eigA =

   -2.8284
    2.0000
    2.8284
   10.0000

Discussion of Results

\( A \) and \( T \) similar: same eigenvalues (not in the same order)
Tridiagonalization took \( n-2 \) multiplications of \( PAP^T \) (plus computation of \( \mathbf{u} \) and \( P = I-2\mathbf{u}\mathbf{u}^T \) in each step).

According to Wolfram Mathworld, “Computing the determinant of such a matrix requires only \( O(7n) \) (as opposed to \( O(\frac{n^3}{3} \))) arithmetic operations) (qtd in Weisstein).”

Householder reflections used in important matrix manipulations, such as QR Factorization.

References

Kreyszig, Erwin, et al. “Chapter 20.9 Tridiagonalization and QR-Factorization.” Advanced Engineering Mathematics: Tenth Edition, 10th ed.,Wiley, Hoboken, NJ, 2011, pp. 888–896.
Weisstein, Eric W. “Tridiagonal Matrix.” From MathWorld–A Wolfram Web Resource. https://mathworld.wolfram.com/TridiagonalMatrix.html
Lee, Che-Rung. “Lecture 4: QR Decomposition.” CS 3331 Numerical Methods. Taiwan, National Tsing Hua University, http://www.cs.nthu.edu.tw/~cherung/teaching/2008cs3331/lec04.pdf. Accessed 26 Apr. 2022.

References continued...

“Orthogonal Matrices: Lecture 8.” Applied Numerical Linear Algebra.
Chalmers University of Technology,
Chalmers University of Technology. http://www.math.chalmers.se/Math/Grundutb/CTH/tma265/1213/Lectures//Lecture8.pdf
Rabinoff, Joseph, and Dan Margalit. Interactive Linear Algebra,
Georgia Institute of Technology, 2017,
textbooks.math.gatech.edu/ila/similarity.html.
Townsend, Alex. “Numerical Linear Algebra and Matrix Factorizations.”
Math 720: Top Ten Algorithms from the 20th Century.
Cornell Department of Mathematics, 2019, Ithica, New York.