LU Decomposition Worked by Hand

Factor a matrix into lower and upper triangular pieces so many linear systems can be solved with cheap back-substitution.

The goal

LU decomposition writes A = L U, with L unit lower-triangular and U upper-triangular. It is Gaussian elimination with the multipliers saved: U is the result of elimination, L records the factors used to zero each column.

3x3 example

Eliminate below the pivot column by column. The multiplier used to clear entry (i,j) becomes L[i,j]; the reduced matrix becomes U. No row of A is discarded - the factors are stored, not thrown away.

import numpy as np
from scipy.linalg import lu
A=np.array([[2.0,1.0,1.0],[4.0,3.0,3.0],[8.0,7.0,9.0]])
P,L,U=lu(A)
print('L=',np.round(L,3)); print('U=',np.round(U,3))
print('check:',np.allclose(P@L@U,A))

Why factor at all

Once you have A = LU, solving A x = b for any right-hand side costs only two triangular solves: forward-substitute L y = b, then back-substitute U x = y - each order n^2. The expensive order n^3 factorization is done once and reused across many b vectors, which is common in time-stepping and design sweeps.

Pivoting

Plain LU fails if a pivot is zero and is unstable if a pivot is tiny. Partial pivoting swaps rows to put the largest available entry on the diagonal, producing PA = LU with a permutation P. This is what production solvers do; it keeps the multipliers bounded and the factorization stable for any nonsingular matrix.