# 7.7: LU Redux

- Page ID
- 1977

\( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)

\( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)

Certain matrices are easier to work with than others. In this section, we will see how to write any square matrix \(M\) as the product of two simpler matrices. We will write $$M=LU\, ,$$ where:

1. \(L\) is \(\textit{lower triangular}\). This means that all entries above the main diagonal are zero. In notation, \(L=(l^{i}_{j})\) with \(l^{i}_{j}=0\) for all \(j>i\).

\[L=\begin{pmatrix}

l^{1}_{1} & 0 & 0 & \cdots \\

l^{2}_{1} & l^{2}_{2} & 0 & \cdots \\

l^{3}_{1} & l^{3}_{2} & l^{3}_{3} & \cdots \\

\vdots & \vdots & \vdots & \ddots \\

\end{pmatrix}

\]

2. \(U\) is \(\textit{upper triangular}\). This means that all entries below the main diagonal are zero. In notation, \(U=(u^{i}_{j})\) with \(u^{i}_{j}=0\) for all \(j<i\).

\[U=\begin{pmatrix}

u^{1}_{1} & u^{1}_{2} & u^{1}_{3} & \cdots \\

0 & u^{2}_{2} & u^{2}_{3} & \cdots \\

0 & 0 & u^{3}_{3} & \cdots \\

\vdots & \vdots & \vdots & \ddots \\

\end{pmatrix}

\]

\(M=LU\) is called an \(\textit{\(LU\) decomposition}\) of \(M\).

This is a useful trick for computational reasons; it is much easier to compute the inverse of an upper or lower triangular matrix than general matrices. Since inverses are useful for solving linear systems, this makes solving any linear system associated to the matrix much faster as well. The determinant---a very important quantity associated with any square matrix---is very easy to compute for triangular matrices.

Example 87

Linear systems associated to upper triangular matrices are very easy to solve by back substitution.

$$\left(\begin{array}{rr|r}

a & b & 1 \\

0 & c & e \\

\end{array}\right) \Rightarrow \ y=\frac{e}{c}\, , \quad x=\frac{1}{a}\left(1-\frac{be}{c}\right)$$

$$\left(\begin{array}{rrr|r}

1 & 0 & 0 & d \\

a & 1 & 0 & e \\

b & c & 1 & f \\

\end{array}\right) \Rightarrow x=d\, , \qquad y=e-ad\, , \qquad z=f-bd-c(e-ad)$$

For lower triangular matrices, \(\textit{back}\) substitution gives a quick solution; for upper triangular matrices, \(\textit{forward}\) substitution gives the solution.

## 7.7.1 Using \(LU\) Decomposition to Solve Linear Systems

Suppose we have \(M=LU\) and want to solve the system

\[ MX=LUX=V.\]

- Step 1: Set \(W=\begin{pmatrix}u\\v\\w\end{pmatrix}=UX\).
- Step 2: Solve the system \(LW=V\). This should be simple by forward substitution since \(L\) is lower triangular. Suppose the solution to \(LW=V\) is \(W_{0}\).
- Step 3: Now solve the system \(UX=W_{0}\). This should be easy by backward substitution, since \(U\) is upper triangular. The solution to this system is the solution to the original system.

We can think of this as using the matrix \(L\) to perform row operations on the matrix \(U\) in order to solve the system; this idea also appears in the study of determinants.

Example 88

Consider the linear system:

$$\begin{array}{r} 6x + 18y + 3z = 3 \\ 2x + 12y + z = 19 \\ 4x + 15y + 3z = 0\end{array}$$

An \(LU\) decomposition for the associated matrix \(M\) is:

\[

\begin{pmatrix}

6 & 18 & 3 \\

2 & 12 & 1 \\

4 & 15 & 3

\end{pmatrix} =

\begin{pmatrix}

3 & 0 & 0 \\

1 & 6 & 0 \\

2 & 3 & 1

\end{pmatrix}

\begin{pmatrix}

2 & 6 & 1 \\

0 & 1 & 0 \\

0 & 0 & 1

\end{pmatrix}.

\]

Step 1: Set \(W=\begin{pmatrix}u\\v\\w\end{pmatrix}=UX\).

Step 2: Solve the system \(LW=V\):

\[

\begin{pmatrix}

3 & 0 & 0 \\

1 & 6 & 0 \\

2 & 3 & 1

\end{pmatrix}

\begin{pmatrix}u\\v\\w\end{pmatrix} =

\begin{pmatrix}3\\19\\0\end{pmatrix}

\]

By substitution, we get \(u=1\), \(v=3\), and \(w=-11\). Then \[W_{0}=\begin{pmatrix}1\\3\\-11\end{pmatrix}\]

Step 3: Solve the system \(UX=W_{0}\).

\[

\begin{pmatrix}

2 & 6 & 1 \\

0 & 1 & 0 \\

0 & 0 & 1

\end{pmatrix}

\begin{pmatrix}x\\y\\z\end{pmatrix} =

\begin{pmatrix}1\\3\\-11\end{pmatrix}

\]

Back substitution gives \(z=-11, y=3\), and \(x=-3\).

Then \(X=\begin{pmatrix}-3\\3\\-11\end{pmatrix}\), and we're done.

## 7.7.2 Finding an \(LU\) Decomposition.

In Section 2.3, Gaussian elimination was used to find \(LU\) matrix decompositions. These ideas are presented here again as review. For any given matrix, there are actually many different \(LU\) decompositions. However, there is a unique \(LU\) decomposition in which the \(L\) matrix has ones on the diagonal. In that case \(L\) is called a \(\textit{lower unit triangular matrix}\).

To find the \(LU\) decomposition, we'll create two sequences of matrices \(L_{1}, L_{2},\ldots\) and \(U_{1}, U_{2}, \ldots\) such that at each step, \(L_{i}U_{i}=M\). Each of the \(L_{i}\) will be lower triangular, but only the last \(U_{i}\) will be upper triangular. The main trick for this calculation is captured by the following example:

Example 89: An Elementary Matrix

Consider $$E=\begin{pmatrix}1&0\\\lambda&1\end{pmatrix}\, ,\qquad M=\begin{pmatrix}a&b&c&\cdots\\d&e&f&\cdots\end{pmatrix}\, .$$

Lets compute \(EM\)

$$EM=\begin{pmatrix}a&b&c&\cdots\\d+\lambda a&e+\lambda b&f+\lambda c&\cdots\end{pmatrix}\, .$$

Something neat happened here: multiplying \(M\) by \(E\) performed the row operation \(R_{2}\to R_{2}+\lambda R_{1}\) on \(M\).

Another interesting fact:

$$E^{-1}:=\begin{pmatrix}1&0\\-\lambda&1\end{pmatrix}$$

obeys (check this yourself...)

$$

E^{-1} E = 1\, .

$$

Hence \(M=E^{-1} E M\) or, writing this out

$$

\begin{pmatrix}a&b&c&\cdots\\d&e&f&\cdots\end{pmatrix}=\begin{pmatrix}1&0\\-\lambda&1\end{pmatrix} \begin{pmatrix}a&b&c&\cdots\\d+\lambda a&e+\lambda b&f+\lambda c&\cdots\end{pmatrix}\, .

$$

Here the matrix on the left is lower triangular, while the matrix on the right has had a row operation performed on it.

We would like to use the first row of \(M\) to zero out the first entry of every row below it. For our running example, $$M=\begin{pmatrix}

6 & 18 & 3 \\

2 & 12 & 1 \\

4 & 15 & 3

\end{pmatrix}\, ,$$ so we would like to perform the row operations $$R_{2}\to R_{2} -\frac{1}{3}R_{1} \mbox{ and } R_{3}\to R_{3}-\frac{2}{3}R_{1}\, .$$

If we perform these row operations on \(M\) to produce

$$U_1=\begin{pmatrix}

6 & 18 & 3 \\

0 & 6 & 0 \\

0 & 3 & 1

\end{pmatrix}\, ,$$

we need to multiply this on the left by a lower triangular matrix \(L_{1}\) so that the product \(L_{1}U_{1}=M\) still.

The above example shows how to do this: Set \(L_{1}\) to be the lower triangular matrix whose first column is filled with \(\textit{minus}\) the constants used to zero out the first column of \(M\). Then $$L_{1} = \begin{pmatrix}

1 & 0 & 0 \\

\frac{1}{3} & 1 & 0 \\

\frac{2}{3} & 0 & 1

\end{pmatrix}\, .$$

By construction \(L_{1}U_{1}=M\), but you should compute this yourself as a double check.

Now repeat the process by zeroing the second column of \(U_{1}\) below the diagonal using the second row of \(U_{1}\) using the row operation \(R_{3}\to R_{3}-\frac{1}{2}R_{2}\) to produce

\[U_{2}=\begin{pmatrix}6&18&3\\0&6&0\\0&0&1\end{pmatrix}\, .\]

The matrix that undoes this row operation is obtained in the same way we found \(L_{1}\) above and is:

$$

\begin{pmatrix}

1&0&0\\

0&1&0\\

0&\frac{1}{2}& 0

\end{pmatrix}\, .

$$

Thus our answer for \(L_{2}\) is the product of this matrix with \(L_{1}\), namely

$$

L_{2}=

\begin{pmatrix}

1 & 0 & 0 \\

\frac{1}{3} & 1 & 0 \\

\frac{2}{3} & 0 & 1

\end{pmatrix}\begin{pmatrix}

1&0&0\\

0&1&0\\

0&\frac{1}{2}& 0

\end{pmatrix}

=\begin{pmatrix}

1 & 0 & 0 \\

\frac{1}{3} & 1 & 0 \\

\frac{2}{3} & \frac{1}{2} & 1

\end{pmatrix}\, .

$$

Notice that it is lower triangular because

$$\textit{The product of lower triangular matrices is always lower triangular!}$$

Moreover it is obtained by recording minus the constants used for all our row operations in the appropriate columns (this always works this way). Moreover, \(U_{2}\) is upper triangular and \(M=L_{2}U_{2}\), we are done! Putting this all together we have

$$M=\begin{pmatrix}

6 & 18 & 3 \\

2 & 12 & 1 \\

4 & 15 & 3

\end{pmatrix}= \begin{pmatrix}

1 & 0 & 0 \\

\frac{1}{3} & 1 & 0 \\

\frac{2}{3} & \frac{1}{2} & 1

\end{pmatrix}\begin{pmatrix}

6 & 18 & 3 \\

0 & 6 & 0 \\

0 & 0 & 1

\end{pmatrix}\, .$$

If the matrix you're working with has more than three rows, just continue this process by zeroing out the next column below the diagonal, and repeat until there's nothing left to do.

The fractions in the \(L\) matrix are admittedly ugly. For two matrices \(LU\), we can multiply one entire column of \(L\) by a constant \(\lambda\) and divide the corresponding row of \(U\) by the same constant without changing the product of the two matrices. Then:

For matrices that are not square, \(LU\) decomposition still makes sense. Given an \(m\times n\) matrix \(M\), for example we could write \(M=LU\) with \(L\) a square lower unit triangular matrix, and \(U\) a rectangular matrix. Then \(L\) will be an \(m\times m\) matrix, and \(U\) will be an \(m\times n\) matrix (of the same shape as \(M\)). From here, the process is exactly the same as for a square matrix. We create a sequence of matrices \(L_{i}\) and \(U_{i}\) that is eventually the \(LU\) decomposition. Again, we start with \(L_{0}=I\) and \(U_{0}=M\).

Example 90

Let's find the \(LU\) decomposition of

\(M=U_{0}=\begin{pmatrix}

-2 & 1 & 3 \\

-4 & 4 & 1

\end{pmatrix}\).

Since \(M\) is a \(2\times 3\) matrix, our decomposition will consist of a \(2\times 2\) matrix and a \(2\times 3\) matrix. Then we start with \(L_{0}=I_{2}=\begin{pmatrix}

1 & 0 \\

0 & 1

\end{pmatrix}\).

The next step is to zero-out the first column of \(M\) below the diagonal. There is only one row to cancel, then, and it can be removed by subtracting \(2\) times the first row of \(M\) to the second row of \(M\). Then:

\[

L_1=\begin{pmatrix}

1 & 0 \\

2 & 1

\end{pmatrix}, \qquad

U_1 = \begin{pmatrix}

-2 & 1 & 3 \\

0 & 2 & -5

\end{pmatrix}

\]

Since \(U_{1}\) is upper triangular, we're done. With a larger matrix, we would just continue the process.

## 7.7.3 Block \(LDU\) Decomposition

Let \(M\) be a square block matrix with square blocks \(X,Y,Z,W\) such that \(X^{-1}\) exists. Then \(M\) can be decomposed as a block \(LDU\) decomposition, where \(D\) is block diagonal, as follows:

\[

M=\begin{pmatrix}

X & Y \\

Z & W

\end{pmatrix}

\]

Then:

$$M=\left(\begin{array}{cc}

I & 0 \\

ZX^{-1} & I

\end{array}\right)

\left(\begin{array}{cc}

X & 0 \\

0 & W-ZX^{-1}Y

\end{array}\right)\left(\begin{array}{cc}

I & X^{-1}Y \\

0 & I

\end{array}\right)\, .$$

This can be checked explicitly simply by block-multiplying these three matrices.

Example 91

For a \(2\times 2\) matrix, we can regard each entry as a \(1\times1\) block.

\[

\begin{pmatrix}

1 & 2 \\

3 & 4

\end{pmatrix}=

\begin{pmatrix}

1 & 0 \\

3 & 1

\end{pmatrix}

\begin{pmatrix}

1 & 0 \\

0 & -2

\end{pmatrix}

\begin{pmatrix}

1 & 2 \\

0 & 1

\end{pmatrix}

\]

By multiplying the diagonal matrix by the upper triangular matrix, we get the standard \(LU\) decomposition of the matrix.