Every essay — page 11
Iterating, instead of factorising
At scale nobody factorises, and the reason is not accuracy — it is that the factors of a sparse matrix are not sparse. What replaces elimination is a sequence of approximations, and the question changes shape: not what the residual of the factorisation is, but how fast the error falls and whether the thing you can measure tracks the thing you cannot. On the model problem every rate is known in closed form before anything runs.
Aggregating what the matrix calls strong
The depth phase measured every method it had on the 45°-rotated anisotropic operator — 0.784, 0.883, 0.844 — and diagnosed the failure as being in the discretisation rather than in the hierarchy.
One sequence and two recurrences
CGLS and LSQR compute the same iterates — the minimiser over a space is unique, so there is nothing to choose between them in the algebra.
Exact along one axis
The tuned diffusion makes the answer exact at every node, and in two dimensions it holds at exactly one flow angle.
The direction the diffusion does not go
Streamline diffusion adds τbbᵀ, a rank-one tensor that annihilates every direction across the flow.
An operator with no entries
At the sizes where linear algebra is expensive the matrix does not exist.
The zero that means it is finished
Every Krylov method ends by dividing by a number the previous step produced, and when that number is zero the recurrence stops.
The same zero, and nothing was found
Change the recurrence by two lines and the divisor stops being a norm.
The division that cannot be done
Conjugate gradients divides by pᵀAp at every step, and on a matrix that is not positive definite that number can be zero or negative.
The residual the method reports
Conjugate gradients prints a relative residual of 6.9·10⁻²¹.
The number that is re-derived
GMRES prints a residual it never computes from its answer either.
An iterate that must be made smaller
Applying a Kronecker-sum operator to a low-rank iterate multiplies its ranks by d and adding two of them adds their ranks, so a solver in a compressed format cannot keep what it produces.
A Krylov space for a problem that is not linear
A quadratic eigenvalue problem has no matrix to build a Krylov space out of.
A different equation on every grid
Upwinding is the exact discretisation of a convection–diffusion problem with diffusion ε + h/2, entry for entry, at a relative difference of between 0 and 1.26·10⁻¹⁶ on every mesh from 15 points to 511.
A parameter that is also a price
ξ = coth(Pe) − 1/Pe is the fraction of h/2 that makes a boundary-layer solution exact at every node.
A run that is over at step five
A conjugate gradient whose every iterate is cut to a rank budget reaches the floor that budget allows at step 5, 36, 42 or 59, and then does nothing for the rest of the run.
A smoother that stops being one
Weighted Jacobi's smoothing factor on the convection–diffusion operator is a function of the cell Péclet number and nothing else — identical to eight digits at five grid sizes at matched Pe.
How much direction there was to lose
At 45° the nine-point stencil hands smoothed aggregation the same wrong hierarchy at every anisotropy — six strong neighbours per interior point, 121 aggregates, the identical partition from ε = 10⁻⁴ to 0.099.
The switch does not know which side is better
The strength threshold moves the coarsening from full to semi at θ = ε exactly, at every anisotropy.
A proof that does not ask how large the matrix is
Proving a Hessian indefinite costs three matrix–vector products when the negative eigenvalue is 3 and nine to eleven when it is a thousandth, and that pair of numbers barely moves across a fourfold range in n.
A walk needs a length
The gap between the two residuals grows as the square root of something, and a square root needs a length.
The answer that arrives when the space runs out
A second-order Krylov recurrence holds vectors of length n for a problem with 2n eigenvalues, so it is exact at n steps where the linearised route needs 2n.
The certificate that arrives soonest is worth least
The more negative a Hessian's smallest eigenvalue, the sooner conjugate gradients meets a direction of negative curvature — and the less of the exact trust-region decrease that direction turns out to be worth.
Sparsity, and what elimination costs
Eliminating a variable couples everything it touched to everything else it touched, and every one of those couplings is an entry that was zero in the matrix and is not zero in its factor. On the same matrix one elimination order gives a factor of a thousand entries and another gives ten thousand — the two factorisations equally accurate, one of them fitting in memory. Nothing numerical chooses between them.