The number that is re-derived

The essay on what a method reports ends with a nine-and-a-half-order gap between a printed residual and a real one, and the natural conclusion to draw from it is the wrong one.

The wrong conclusion is that iterative methods print residuals that cannot be trusted. The right one is narrower and far more useful, and getting to it requires a method in the same family that does not have the problem.

GMRES is that method. It also prints a residual it never computes from its answer. Its number is honest.

Where GMRES’s number comes from #

GMRES builds an Arnoldi basis and solves a small least-squares problem in it. The residual it reports is not updated by a recurrence at all: at each step, the Givens rotations that triangularise the Hessenberg matrix carry the right-hand side along with them, and the last entry of the rotated right-hand side is the residual norm of the projected problem. Exactly, as a consequence of the rotations being orthogonal.

So the quantity is re-derived at every step from the whole of what the method currently holds, rather than updated from its own previous value. That is the difference from a conjugate gradient residual, and it is the only difference that matters here.

Three matrices chosen to break it #

The comparison is only worth having if the matrices are hostile, so all three are.

The diagonal one with an eigenvalue at 10⁻¹⁴. The same matrix the recurrence essay is built on, where the conjugate gradient recurrence goes 7.3·10¹⁰ wrong. Worst reported-to-actual factor: 1.80.

A bidiagonal matrix with a superdiagonal thirty times its diagonal. Strongly non-normal, the kind of operator whose eigenvalues predict nothing about its behaviour. Worst factor: 2.86.

A single Jordan block. The most degenerate matrix in the subject: one eigenvalue, one eigenvector, and a Krylov space that cannot resolve it. Worst factor: 2.00.

Not one of the three produces a discrepancy worth a sentence. The largest is a factor of under three, at one step out of forty, and the two curves are visually one curve on all three matrices.

While the basis is gone #

Here is the part that makes this a measurement rather than a reassurance.

The Arnoldi basis is built by modified Gram–Schmidt with no reorthogonalisation, which is what GMRES does. On two of the three matrices, ‖VᵀV − I‖ after forty steps is 1.41.

That is not a loss of orthogonality, it is the absence of it: a value of order one means the basis vectors are no longer even approximately mutually orthogonal, and the projection GMRES is computing is a projection onto something other than the space it thinks it is projecting onto.

And the residual it derives from that projection is still within a factor of two of the truth.

The reason is a theorem rather than luck — modified-Gram–Schmidt GMRES is backward stable, and the loss of orthogonality is exactly compensated by what it does to the projected problem — but the measurement is what this collection is for, and the measurement says: the basis is destroyed and the number is fine.

And it is not forty steps being too few #

A basis at ‖VᵀV − I‖ = 1.41 has stopped being a basis, and the reported residual is still within a factor of two. That invites one obvious objection: forty steps may simply not be long enough for the trouble to arrive, and the measurement is of a run that has not yet gone wrong.

Taking each matrix to n steps at three sizes, worst reported-to-actual factor with the orthogonality loss beside it:

tiny gives 1.80 at n = 40 with a loss of 1.41, 4.47 at n = 100 with 1.92, and 1.48 at n = 200 with 2.51. bidiagonal gives 2.86 with 0.00, 3.04 with 2.00, and 3.67 with 2.46. jordan gives 2.00 with 1.41, 1.90 with 0.00, and 2.58 with 1.41.

The worst factor over all nine runs is 4.47, and it has no trend in the size at all — on tiny the largest is at a hundred steps and the smallest at two hundred.

The orthogonality loss does have a trend. It climbs monotonically on two of the three matrices: 1.41, 1.92, 2.51 on tiny, and 0.00, 2.00, 2.46 on bidiagonal. So at two hundred steps the basis is measurably further from orthogonal than it was at forty, and the number derived from it is no further from the truth.

Which is the whole of the pair, stated as a growth rate #

That closes the comparison in the form that makes it a law rather than an anecdote.

A carried quantity’s error grows with the number of steps it has been carried through. The conjugate gradient recurrence’s gap is a random walk — the essay before this one fits it at 0.507 — so it climbs like the square root of the run, without limit, and a longer run is strictly worse.

A re-derived quantity’s error does not grow with anything. It has no history to accumulate: not the roundings of the previous steps, and — the measurement here — not even the history of the basis it is derived from having fallen apart. Nine runs, five times the step count, a basis degrading monotonically, and one bounded factor.

So the two methods do not differ by how carefully each is implemented, or by how well conditioned each problem is, or by how orthogonal each basis stays. They differ in whether the printed number is a function of the current state or a function of every state so far, and that difference is a growth rate against the run length: √k for one and nothing for the other.

That is also the practical rule it leaves. A quantity a method reports should be recomputed from what the method currently holds, wherever recomputing it is affordable — and the interesting half is how often that is. GMRES’s costs nothing: it falls out of rotations the method applies anyway. Conjugate gradients’ costs one matrix–vector product, which is one in a hundred and fifty over a run and buys a factor of 7·10¹⁰. Neither of those is a difficult trade, and the reason the defect survives is not that the repair is expensive.

What the pair actually says #

Two methods. The same family, the same Krylov space, the same matrix, the same arithmetic. One of them reports a residual nine and a half orders out and the other reports one within a factor of two.

The difference is not accuracy per step. It is not the conditioning, which is identical. It is not the orthogonality, which is worse in the method that gets it right. It is:

one method carries a vector forward and the other re-derives a number.

A carried quantity accumulates every rounding that has ever entered it, and nothing downstream will ever compare it against anything. A re-derived quantity has no history: whatever was wrong with the last step’s version of it is not in this step’s, because this step’s was computed from scratch.

That is the sentence this field is organised around, and it is worth having it in the form where the two sides are two implementations of the same idea rather than two different subjects.

The price GMRES pays for it #

An honest comparison has to say what the re-derivation costs, because it is not free and the reason conjugate gradients uses a recurrence is a good one.

GMRES keeps the whole basis. Storage grows with the step number and so does the work per step: at step k it orthogonalises against k vectors. That is why it is restarted in practice, and restarting brings its own essay’s worth of problems.

Conjugate gradients keeps three vectors and does a fixed amount of work per step, for as many steps as anybody wants. On a symmetric matrix that is an enormous prize, and it is bought exactly by the short recurrences — including the one for the residual.

So the two methods are not a good one and a bad one. They are two points on a trade this collection has drawn before: a short recurrence is cheap and forgets nothing, and a long one is expensive and starts fresh. The residual gap is one more thing on the bill for the short recurrence, and it is usually the smallest item on it.

What “re-derived” means, precisely #

The word is doing a lot of work in this essay, so it is worth pinning down what distinguishes the two cases, because “computed from scratch” is not quite it.

GMRES’s residual at step k is not computed from scratch either. The Givens rotations are applied incrementally — the rotations from earlier steps are stored and reapplied — so there is a recurrence in the implementation. What is different is what the recurrence is on.

The conjugate gradient recurrence is on the residual vector itself: the quantity being reported is the object being updated, so an error in it is an error in the report, permanently.

The GMRES recurrence is on the factorisation of the projected problem. The reported number is the last entry of a rotated right-hand side, and it is a consequence of a triangularisation whose own errors are bounded by the orthogonality of the rotations — which is exact to a rounding, whatever has happened to the Arnoldi basis. An error in the basis changes which problem is being projected onto, and the reported residual is then the correct residual of that slightly different problem — which is precisely a backward-error statement, and precisely why the number stays honest.

So the distinction is not “recurrence or not”. It is whether the reported quantity is the state of a recurrence, or a function of a state whose own errors are backward errors. The first accumulates into the report. The second accumulates into the problem, where this collection has a whole vocabulary for it.

The cost of asking, once #

There is one number this pair of essays has not priced and it is the one a reader would act on, so it is worth putting down.

Recomputing ‖b − Ax‖ from the returned vector is one matrix–vector product. On a run of a hundred and fifty iterations that is two thirds of one per cent of the work, paid once, at the end, after the iteration has stopped. It does not change the answer, it does not change the search directions, and it cannot fail.

What it buys is that the number printed beside an answer is a number about that answer. On the conjugate gradient runs measured here the difference between the two is a factor of 7.3·10¹⁰; on almost every ordinary run it is a factor of one, and the product was wasted. That asymmetry is the whole argument: a check that is free when it is unnecessary and decisive when it is not is a check worth making unconditionally, which is the same reasoning the collection’s own rule about drawing a residual beside a factorisation rests on.

And a third quantity, which is neither #

There is a number GMRES does not report and could not, and it is worth naming so the essay does not read as an endorsement.

GMRES’s residual is honest about the vector it would return. It says nothing about the error, and the collection’s standing warning applies here exactly as it does everywhere else: a relative residual of 10⁻¹⁵ on the Jordan block is a statement about ‖b − Ax‖ and not about ‖x − x*‖, and on a matrix with a condition number of 10¹⁰ the second can be ten orders larger than the first.

So the ledger for a GMRES run has three quantities in it, and the essays in this field have now separated all three:

• what the method printed — honest here, and not in the conjugate gradient case;
• the residual of the answer — which is what the printed number should equal, and does;
• the error of the answer — which neither quantity measures and which the condition number amplifies from the first.

Which methods are on which side #

The pair in this essay is two points, and it is worth saying what the rest of the family looks like, because the classification is usable rather than anecdotal.

Methods whose reported residual is carried. Conjugate gradients, the whole family of short recurrences built on it, and every method whose implementation updates a residual vector because it cannot afford a second product: BiCG, BiCGStab, CGS, and conjugate gradients on the normal equations. All of them have a residual that is the state of a recurrence, all of them are subject to the mechanism the recurrence essay measures, and the size of the effect is set by how far the iterates travel rather than by anything about the method.

Methods whose reported residual is re-derived. GMRES and its relatives, MINRES, and anything else that reports the residual of a projected problem obtained by an orthogonal factorisation. Their number is a backward-error statement about a nearby problem, and the loss of orthogonality moves which nearby problem rather than corrupting the number.

The classification is not about symmetry, it is not about storage, and it is not about which method is better. It is about a single implementation decision that is usually made for cost reasons and is not usually written down as an accuracy decision at all.

And the useful consequence is the one a reader can act on: if a solver’s reported residual comes from a short recurrence, the number it prints is a claim about the run and not about the answer. One extra product at the end converts the claim; nothing else needs to change.

What a reader should take from the pair #

The two essays together are a claim about implementations rather than about methods, and it is worth separating what a reader can act on from what is merely interesting.

What is interesting is the mechanism: a residual carried in a vector accumulates every rounding that enters it, and a residual re-derived from an orthogonal factorisation does not. That is a fact about arithmetic and it does not depend on anybody’s code.

What is actionable is narrower. If a solver reports a residual and the method behind it uses a short recurrence, the number it prints is a claim about the run rather than about the answer, and one extra product converts it into the second. If the method holds a full basis and derives the number from a projection, the printed number is already a statement about the vector that will be returned, and there is nothing to do.

Neither of those requires knowing the theory. Both require knowing which method is underneath, which is the one thing a caller usually does know and the one thing a residual printed on its own does not say. The collection’s rule about drawing a factorisation with its residual attached is the same idea one level down: a number and the thing it describes should not be separable by a reader skimming, and here they are separated by an implementation decision nobody wrote down as one.

The refusal #

The claim under test is the over-general one: that a Krylov method’s reported residual cannot be trusted.

It is fed the conjugate gradient run on the same matrix, where the factor between the reported number and the answer’s is 7.3·10¹⁰. Fed that, the assertion that GMRES’s worst factor is under five is false about the wrong method, and the paired assertion — that on the same matrix the recurrence is a million times further out than the projection ever is — passes.

The two assertions together are the essay. Neither of them alone says anything: “GMRES is accurate” is a reassurance, and “conjugate gradients is not” is a scare. The measurement is the ratio between them on one matrix, and it is 7.3·10¹⁰ against 1.80.

A worst factor of 1.80 on a basis that has ceased to be orthogonal is the essay’s whole claim in one frame. The second matrix says it is not a property of that particular defect:

The third matrix is the one that breaks the pattern a reader will have started forming.

The residual is least honest where the basis is most orthogonal, which is the opposite of what the essay’s setup invites a reader to expect. The two matrices that destroy orthogonality entirely give factors of 1.80 and 2.00; the one that preserves it to nine digits gives 2.86. Three matrices is not a correlation, but it is enough to rule out the reading that the Givens residual degrades because the basis does — on this evidence it does not degrade with the basis at all.

And 1.41 is √2, twice. That is not a coincidence and it is not a coarse reading: ‖VᵀV − I‖ reaches √2 when two of the basis vectors have become parallel, so it is the saturation value rather than a large number. Both of the first two matrices are at it. The quantity has stopped being able to report how much worse things are getting, which is worth knowing before comparing two runs by it.

So the essay’s claim survives its own hardest case and survives it for a reason the three matrices make visible. The Givens residual is re-derived from the whole basis at every step, so what it inherits is the basis’s span rather than the basis’s orthogonality — and a span can be right while the vectors describing it are not. A carried residual inherits the arithmetic that produced it and has no such recourse.