Section 1.8 \(LU\) Decomposition
The approach we used in
Theorem 1.7.12 also gives us a useful way to
factor matrices. As with polynomials, the goal in factoring a matrix
\(A\) is to write it as the product of simpler matrices. The factorization we'll look at in this section is the
\(LU\) decomposition, which we demonstrate first by an example.
Example 1.8.1. \(LU\) Factorization.
Let
\begin{equation*}
A = \mqty[3 & 2 & -1 \\ -2 & 4 & 4 \\ 3 & 1 & 2]\text{.}
\end{equation*}
Now consider the problem of reducing \(A\) to an echelon form. First, we can take \(\frac{2}{3}\) times row 1 and add it to row 2 and \(-1\) times row 1 and add it to row 3 to obtain
\begin{equation*}
\mqty[3 & 2 & -1 \\ 0 & \frac{16}{3} & \frac{10}{3} \\ 0 & -1 & 3].
\end{equation*}
A final operation, adding \(\frac{3}{16}\) times row 2 to row 3, gives us the echelon form
\begin{equation*}
U = \mqty[3 & 2 & -1 \\ 0 & \frac{16}{3} & \frac{10}{3} \\ 0 & 0 & \frac{29}{8}]
\end{equation*}
Now, we also saw in the last section that row operations correspond to multiplication by elementary matrices. In this case, \(U = E_{23}E_{13}E_{12}A\) where
\begin{equation*}
E_{12} = \mqty[1 & 0 & 0 \\ \frac{2}{3}& 1 & 0 \\ 0 & 0 & 1], E_{13} = \mqty[1 & 0 & 0 \\ 0 & 1 & 0 \\ -1 & 0 & 1]\text{ and }E_{23} = \mqty[1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & \frac{3}{16} & 1].
\end{equation*}
If we now solve the equation \(U = E_{23}E_{13}E_{12}A\) for \(A\text{,}\) we get the useful product
\begin{equation*}
A = E_{12}^{-1}E_{13}^{-1}E_{23}^{-1}U.
\end{equation*}
It turns out to be very easy to find \(E_{12}^{-1}E_{13}^{-1}E_{23}^{-1}\) as these are elementary matrices, and we get
\begin{equation*}
E_{12}^{-1}E_{13}^{-1}E_{23}^{-1} = \mqty[1 & 0 & 0 \\ -\frac{2}{3} & 1 & 0 \\ 1 & -\frac{3}{16} & 1]
\end{equation*}
This matrix is lower triangular, and we denote it by \(L\text{.}\) Thus we have the factorization \(A = LU\text{.}\)
We can verify the factorization in
Example 1.8.1 using Octave as below:
The Octave command
lu() can also be used to find matrices
\(L\) and
\(U\) satisfying
\(A = LU\text{:}\)
The matrix factorization \(A = LU\) demonstrated above can always be computes (reordering the rows of \(A\) if necessary). In this factorization, \(L\) is lower triangular with \(1\)s along the diagonal and \(U\) is an echelon form of \(A\) (and therefore upper triangular if \(A\) is square). This form has several advantages. First, if \(A\) is square then \(\det(A) = \det(L)\det(U) = \det(U)\text{,}\) and the determinant of \(U\) is very easy to find. Second, this form lets us solve systems more quickly.
Example 1.8.2. Solving a System with \(LU\) factorization.
Let
\begin{equation*}
A = \mqty[-1 & -7 & -3 \\ 2 & 15 & 6 \\ 1 & 3 & 2]\text{ and }\vb{b} = \mqty[1 \\ -5 \\ 11].
\end{equation*}
Find an \(LU\) decomposition of \(A\) and use it to solve \(A\vb{x} = \vb{b}\text{.}\)
Solution.
We can find \(L\) and \(U\) by reducing \(A\) to an echelon form. One possibility is
\begin{equation*}
L = \mqty[1 & 0 & 0 \\ -2 & 1 & 0 \\ -1 & -4 & 1]\text{ and }U = \mqty[-1 & -7 & -3 \\ 0 & 1 & 0 \\ 0 & 0 & -1].
\end{equation*}
Now we can solve \(A\vb{x} = \vb{b}\) in two steps:
Solve \(L\vb{y} = \vb{b}\) for \(\vb{y}\text{.}\)
Solve \(U\vb{x} = \vb{y}\) for \(\vb{x}\text{.}\)
Solving the first equation by reducing \(\smqty[L & \vb{b}]\) produces
\begin{equation*}
\vb{y} = \mqty[1 \\ -3 \\ 0]\text{.}
\end{equation*}
Then, solving the second equation by reducing \(\smqty[U & \vb{y}]\) produces
\begin{equation*}
\vb{x} = \mqty[20\\-3\\0].
\end{equation*}
The process in
Example 1.8.2 produces the same exact answer that row reduction would if applied to
\(\smqty[A & \vb{b}]\text{,}\) or if
\(A^{-1}\) were found and multiplied to
\(\vb{b}\text{.}\) However, computing and using the
\(LU\) factorization can be more stable numerically than either Gaussian elimination or computing the inverse.
