A Tour of Matrix Decompositions — When to Reach for LU, QR, Cholesky, or SVD

A decomposition rewrites one awkward transformation as a product of parts you can actually work with: solve it (LU), solve it at half price (Cholesky), orthogonalize and fit (QR), or survive a broken matrix (SVD). Starting from why nobody computes an inverse, this piece builds up to a single table that tells you which one to reach for, plus the numpy that goes with it.

Solving without ever forming an inverse

Ask anyone how to solve $Ax = b$ and they will write $x = A^{-1}b$. Yet in numerical code, the function that builds an inverse is almost never called. Library documentation says it outright: use solve, not inv.

There are two reasons. Computing a full inverse costs distinctly more than solving the system once. And routing the answer through an inverse mixes in extra error along the way. You want a single vector $x$, but you have built an entire $n \times n$ table first and then multiplied by it — a detour when a straight road exists.

So what do you do instead? You decompose. Rewriting 60 as $2^2 \cdot 3 \cdot 5$ makes its divisors obvious; rewriting a matrix as a product of easy pieces makes everything downstream easy in the same way. LU, Cholesky, QR and SVD are all the same idea in four costumes. What separates them is only this: what you are willing to assume about the matrix, and how cheap or how sturdy you get in return.

The previous article, The Linear Algebra Behind LoRA and RAG, treated a matrix as a machine that deforms space. This one is about taking that deformation apart so you can carry it.

There are really only two kinds of easy piece

The destinations of decomposition are surprisingly few. Triangular matrices, and orthogonal ones.

A triangular matrix has zeros everywhere below (or above) the diagonal. What makes that valuable is that the equations fall like dominoes.

Put in words: the $i$-th answer is fixed by substituting the answers you already have and subtracting. The first row has one unknown, so it resolves immediately. The second row, once you plug in the first answer, also has one unknown left. And so on down. It is called a system of equations, but in practice it is top-to-bottom substitution. Sweeping downward is forward substitution, upward is back substitution, and either costs on the order of $n^2$ — not $n^3$.

An orthogonal matrix has columns that are mutually perpendicular and all of length one. That buys you two things.

Two short claims, which say in words that the transpose is the inverse, and that multiplying by $Q$ never changes a length. The first means you never have to compute an inverse at all. The second means error does not get inflated every time you multiply. That second property is why numerical analysts are so fond of orthogonal matrices: they refuse to amplify what is already wrong.

Perpendicular just means the dot product is zero. Spin the two vectors below and watch the value pass through zero exactly at a right angle. The $Q$ that comes out of a decomposition satisfies that relationship for every pair of its columns at once.

So a matrix decomposition is really this: rewriting a matrix as a product of pieces that solve like dominoes and pieces that do not amplify error. The four decompositions below are just different ways of mixing those two ingredients.

LU — elimination, with the receipts kept

The first one is the elimination you did in school. Double row one, subtract it from row two, repeat until everything below the diagonal is zero. Hand that to a computer and keep the multipliers instead of throwing them away, and you have LU.

$L$ is lower triangular with ones on the diagonal (it holds the multipliers you eliminated with), $U$ is what elimination left above the diagonal, and $P$ is a record of which rows got swapped. Put in words: shuffle the rows suitably and any square matrix can be written as lower triangular times upper triangular.

$P$ earns its place because the divisor at each step — the pivot — breaks everything if it is zero or merely tiny. The fix is to bring the row with the largest absolute value in that column to the top before dividing. That is partial pivoting, and every production LU has it. Most of the ways textbook elimination falls apart in real code trace back to this one spot, and the entire remedy is compressed into that single letter $P$.

Once decomposed, solving is two sweeps: solve $Ly = Pb$ forward, then $Ux = y$ backward. The decomposition costs on the order of $n^3$; the two sweeps cost $n^2$.

This is where the practical payoff lives. When $A$ stays put and only $b$ changes, you decompose once. A hundred right-hand sides, one expensive step. Faster than building an inverse and multiplying a hundred times, and more accurate too. That is all "use solve, not inv" is really saying.

It is worth having a feel for the cost. LU needs roughly $\tfrac{2}{3}n^3$ multiply-divides. Double $n$ and that is eight times the work; ten times $n$ is a thousand times the work. At $n=1000$ it is a few hundred million operations, under a second on a current CPU. At $n=10{,}000$ it is a thousand times that, and only there do you start questioning whether solving a dense matrix at all was the right design. Forward and back substitution, meanwhile, stay at $n^2$ — a mere million operations at $n=1000$. The decomposition is the expensive part; the solving is nearly free, and that gap between the two exponents is exactly why reusing a factorization pays.

LU throws in one more thing: the determinant. $\det A = \pm\prod_i u_{ii}$ — in words, multiply the diagonal of $U$ and let the parity of the row swaps decide the sign. Cofactor expansion by the definition has a number of terms proportional to $n!$ and is hopeless by $n=20$, but through LU the determinant falls out as a byproduct of work you already did. Open up your library's determinant function and this is usually what is inside.