Gaussian elimination

Imagine a stack of linear equations in several unknowns, like a tangle of strings each tying the unknowns together. Gaussian elimination is the patient way of untangling them: you use one equation to wipe a chosen unknown out of all the others, then move to the next, until the last equation involves only one unknown — which you can read off — and then you walk backwards filling in the rest. It is the method you were taught for 'solving by adding and subtracting equations', written down as a careful, repeatable recipe.

Concretely, write the system as A x = b and lay A and b side by side as an augmented array. For column 1 you pick the diagonal entry a_11 as the pivot, and for each row i below you subtract (a_i1 / a_11) times row 1 from row i, which makes every entry under the pivot zero. Now column 1 is cleared except at the top; repeat on column 2 using a_22 as the pivot, and so on. After n-1 such sweeps A has become an upper-triangular matrix U (everything below the diagonal is zero), and the same operations applied to b give a transformed right-hand side. The triangular system is then solved by back substitution. A tiny example: from x + y = 3 and 2x + y = 4, subtract 2 times the first row from the second to get -y = -2, so y = 2, then back-substitute x = 1.

This is the workhorse beneath almost everything in numerical linear algebra; the elimination it performs is exactly what the LU factorization records, so the two are the same computation seen from two angles. It costs about 2 n^3 / 3 floating-point operations for an n-by-n system. The naive version can fail or lose accuracy if a pivot is zero or tiny, which is why in practice it is always done with pivoting (row swaps) for stability.