A calculus of the absurd

22.2.4 Performing LU decomposition

This is all very well, but how do we split $A$ into $LU = A$ (where $L$ is lower-triangular and $U$ is upper-triangular)? Well, we can obtain an upper-triangular matrix by using Gaussian elimination, but what about the lower-triangular matrix? We can represent the series of elementary row operations which we need to perform in order to write $A$ as an upper-triangular matrix as $E_{n}E_{n-1} ... E_{1}$. We also know that these operations have an inverse (the set of matrices which "undo" all the elementary row operations which we have applied)! If we are lucky, we didn’t have to swap any rows, and this inverse is a lower-triangular matrix. If not, then we cannot write $A = LU$ where $L$ is a lower-triangular matrix, but we can write $PA = LU$. Why?

After Gaussian elimination, if we carried out steps represented by the matrix $E$, then we have

\begin{align} A &= IA \\ &= (E^{-1} E) A && \text {As any matrix multiplied by its inverse is $I$} \\ &= (E^{-1}) EA && \text {As matrix multiplication is associative} \‚ \end{align}

Now we have an equation in the form

\begin{align} A = (E^{-1}) EA \end{align}

Unfortunately, $E^{-1}$ is not a lower-triangular matrix (we considered the case where it is above), but what we can do is rearrange some of the rows to make it into one, i.e. we have a matrix $P$ which swaps rows (a “transposition matrix”) of a matrix. In this case,

\begin{align} PA = P(E^{-1})EA \end{align}

That is, we have

\begin{equation} PA = LU \end{equation}

Where $L=PE^{-1}$ and $U=EA$.