LU Decomposition

LU decomposition is a reformulation of Gaussian elimination that makes the factorization reusable. Instead of discarding the row operations after solving one system, we encode them into two triangular matrices $L$ (lower) and $U$ (upper) such that $A = LU$. Solving $Ax = b$ then splits into two triangular systems — $Ly = b$ (forward substitution) followed by $Ux = y$ (back substitution) — each of which costs only $n^2$ operations.

---

Frobenius Matrices

Each round of Gaussian elimination can be represented as multiplication by a Frobenius matrix (also called an elementary lower triangular matrix). A Frobenius matrix is a unit lower triangular matrix that differs from the identity only in one column below the diagonal.

For the $k$-th elimination step, the Frobenius matrix $F^{(k)}$ places the negated multipliers $-m_{ik}$ (for $i > k$) in column $k$:

$F^{(k)} = \begin{bmatrix}1&&&&&\\&\ddots&&&&\\&&1&&&\\&&-m_{k+1,k}&1&&\\&&\vdots&&\ddots&\\&&-m_{n,k}&&&1\end{bmatrix}$

Proof: $A_{k+1} = F^{(k)} A_k$

Let $A_k$ be the matrix at the start of step $k$, with entries $a^{(k)}_{ij}$ and rows $R_1, \dots, R_n$, and let $m_{ik} = a^{(k)}_{ik} / a^{(k)}_{kk}$. We show that left-multiplying by $F^{(k)}$ performs exactly the row operations $R_i \leftarrow R_i - m_{ik} R_k$ for every $i > k$.

The $i$-th row of a product $FA$ is the $i$-th row of $F$ times $A$, which is a combination of the rows of $A$ weighted by the entries of that row of $F$:

$\text{row}_i(FA) = f_{i,1}R_1 + f_{i,2}R_2 + \cdots + f_{i,n}R_n$

Now read the rows of $F^{(k)}$:

• For $i \le k$, row $i$ of $F^{(k)}$ is just row $i$ of the identity, so row $i$ of $F^{(k)}A_k$ is $R_i$. The first $k$ rows are unchanged.
• For $i > k$, row $i$ of $F^{(k)}$ has $-m_{ik}$ in column $k$ and $1$ in column $i$, so row $i$ of $F^{(k)}A_k$ is $R_i - m_{ik}R_k$.

This is the $k$-th round of Gaussian elimination. The column-$k$ entry of each row $i > k$ becomes

$a^{(k)}_{ik} - m_{ik}\,a^{(k)}_{kk} = a^{(k)}_{ik} - \frac{a^{(k)}_{ik}}{a^{(k)}_{kk}}\,a^{(k)}_{kk} = 0$

and the zeros already in columns $1, \dots, k-1$ stay zero, because $R_k$ is also zero in those columns. So $F^{(k)}A_k = A_{k+1}$. $\blacksquare$

Example: $A_2 = F^{(1)}A_1$ for a $3 \times 3$ matrix. With pivot $a_{11}$, the multipliers are $m_{21} = a_{21}/a_{11}$ and $m_{31} = a_{31}/a_{11}$:

$F^{(1)}A_1 = \begin{bmatrix}1&0&0\\-m_{21}&1&0\\-m_{31}&0&1\end{bmatrix}\begin{bmatrix}a_{11}&a_{12}&a_{13}\\a_{21}&a_{22}&a_{23}\\a_{31}&a_{32}&a_{33}\end{bmatrix}$

Taking one row of $F^{(1)}$ at a time:

$\begin{bmatrix}1&0&0\end{bmatrix}A_1 = R_1 \qquad \begin{bmatrix}-m_{21}&1&0\end{bmatrix}A_1 = R_2 - m_{21}R_1 \qquad \begin{bmatrix}-m_{31}&0&1\end{bmatrix}A_1 = R_3 - m_{31}R_1$

Written out entry by entry:

$F^{(1)}A_1 = \begin{bmatrix}a_{11}&a_{12}&a_{13}\\a_{21}-m_{21}a_{11}&a_{22}-m_{21}a_{12}&a_{23}-m_{21}a_{13}\\a_{31}-m_{31}a_{11}&a_{32}-m_{31}a_{12}&a_{33}-m_{31}a_{13}\end{bmatrix}$

The multipliers were chosen so that the first-column entries below the diagonal cancel:

$a_{21}-m_{21}a_{11} = a_{21}-\frac{a_{21}}{a_{11}}\,a_{11} = 0 \qquad a_{31}-m_{31}a_{11} = a_{31}-\frac{a_{31}}{a_{11}}\,a_{11} = 0$

$A_2 = F^{(1)}A_1 = \begin{bmatrix}a_{11}&a_{12}&a_{13}\\0&a^{(2)}_{22}&a^{(2)}_{23}\\0&a^{(2)}_{32}&a^{(2)}_{33}\end{bmatrix}, \qquad a^{(2)}_{ij} = a_{ij} - m_{i1}\,a_{1j}$

Step 2 of the worked example below is this same computation with numbers.

From $A$ to $U$

Applying all $n-1$ Frobenius matrices in turn ($A_1 = A \to A_2 \to \cdots \to A_n = U$) transforms $A$ into $U$:

$F^{(n-1)}\cdots F^{(2)}F^{(1)}A = U$

The lower triangular factor $L$ is then recovered by inverting the Frobenius matrices in order:

$L = \left(F^{(1)}\right)^{-1}\left(F^{(2)}\right)^{-1}\cdots\left(F^{(n-1)}\right)^{-1}$

Note

The inverse of a Frobenius matrix is obtained simply by flipping the signs of its off-diagonal entries. For $F^{(k)}$ with entries $-m_{ik}$ below the diagonal, $\left(F^{(k)}\right)^{-1}$ has entries $+m_{ik}$ in exactly the same positions. This makes assembling $L$ trivial: just collect all multipliers $m_{ik}$ at their natural positions.

---

Worked Example

Problem: Decompose the matrix $A$ below into $L$ and $U$, then use LU decomposition to solve $Ax = b$.

$A=\begin{bmatrix}1&2&1\\1&-2&2\\2&12&-2\end{bmatrix},\quad b=\begin{bmatrix}0\\4\\4\end{bmatrix}$

Step 1 — First Frobenius matrix $F^{(1)}$

Label the original matrix $A_1 = A$. Compute the multipliers for the first elimination column using the pivot $a_{11} = 1$:

$m_{21} = \frac{a_{21}}{a_{11}} = \frac{1}{1} = 1 \qquad m_{31} = \frac{a_{31}}{a_{11}} = \frac{2}{1} = 2$

Build $F^{(1)}$ by placing $-m_{21}$ and $-m_{31}$ in column 1 below the diagonal:

$F^{(1)}=\begin{bmatrix}1&0&0\\-1&1&0\\-2&0&1\end{bmatrix}$

Step 2 — Compute $A_2 = F^{(1)} A_1$

$A_{2}=F^{(1)}A_{1}=\begin{bmatrix}1&0&0\\-1&1&0\\-2&0&1\end{bmatrix}\begin{bmatrix}1&2&1\\1&-2&2\\2&12&-2\end{bmatrix}=\begin{bmatrix}1&2&1\\0&-4&1\\0&8&-4\end{bmatrix}$

The first column below the diagonal is now zero, exactly as in the Gaussian elimination step.

Step 3 — Second Frobenius matrix $F^{(2)}$

Use $A_2$ to compute the multiplier for the second elimination column. The new pivot is $a_{22} = -4$:

$m_{32} = \frac{a_{32}}{a_{22}} = \frac{8}{-4} = -2$

Build $F^{(2)}$ by placing $-m_{32} = 2$ in position $(3,2)$:

$F^{(2)}=\begin{bmatrix}1&0&0\\0&1&0\\0&2&1\end{bmatrix}$

Step 4 — Compute $U = F^{(2)} A_2$

$U=F^{(2)}A_{2}=\begin{bmatrix}1&0&0\\0&1&0\\0&2&1\end{bmatrix}\begin{bmatrix}1&2&1\\0&-4&1\\0&8&-4\end{bmatrix}=\begin{bmatrix}1&2&1\\0&-4&1\\0&0&-2\end{bmatrix}$

This is the upper triangular factor $U$.

Step 5 — Assemble $L$

Rather than explicitly inverting the Frobenius matrices, we reuse the multipliers from each step and place them directly at their natural positions:

$L=\begin{bmatrix}1&0&0\\m_{21}&1&0\\m_{31}&m_{32}&1\end{bmatrix}=\begin{bmatrix}1&0&0\\1&1&0\\2&-2&1\end{bmatrix}$

Verification: $L \cdot U$ should equal $A$:

$\begin{bmatrix}1&0&0\\1&1&0\\2&-2&1\end{bmatrix}\begin{bmatrix}1&2&1\\0&-4&1\\0&0&-2\end{bmatrix}=\begin{bmatrix}1&2&1\\1&-2&2\\2&12&-2\end{bmatrix} = A \checkmark$

---

Solving $Ax = b$ via $LU$

With the factorization $A = LU$, the system $Ax = b$ becomes $LUx = b$. Define $y = Ux$, so that $Ly = b$ and then $Ux = y$.

Step 6 — Forward substitution: solve $Ly = b$

$\begin{bmatrix}1&0&0\\1&1&0\\2&-2&1\end{bmatrix}\begin{bmatrix}y_{1}\\y_{2}\\y_{3}\end{bmatrix}=\begin{bmatrix}0\\4\\4\end{bmatrix}$

Row 1: $y_1 = 0$

Row 2: $y_1 + y_2 = 4 \implies 0 + y_2 = 4 \implies y_2 = 4$

Row 3: $2y_1 - 2y_2 + y_3 = 4 \implies 0 - 8 + y_3 = 4 \implies y_3 = 12$

$y = \begin{bmatrix}0\\4\\12\end{bmatrix}$

Step 7 — Back substitution: solve $Ux = y$

$\begin{bmatrix}1&2&1\\0&-4&1\\0&0&-2\end{bmatrix}\begin{bmatrix}x_{1}\\x_{2}\\x_{3}\end{bmatrix}=\begin{bmatrix}0\\4\\12\end{bmatrix}$

Row 3: $-2x_3 = 12 \implies x_3 = -6$

Row 2: $-4x_2 + x_3 = 4 \implies -4x_2 + (-6) = 4 \implies -4x_2 = 10 \implies x_2 = -2.5$

Row 1: $x_1 + 2x_2 + x_3 = 0 \implies x_1 + 2(-2.5) + (-6) = 0 \implies x_1 = 11$

Step 8 — Solution

$\boxed{x_1 = 11,\quad x_2 = -2.5,\quad x_3 = -6}$

Note

The power of LU decomposition is that $L$ and $U$ only need to be computed once. If we need to solve $Ax = b$ for many different right-hand side vectors $b$ (with the same $A$), we simply repeat Steps 6 and 7 with the new $b$ each time. This is far cheaper than re-running Gaussian elimination from scratch for every new $b$.