Doing it twice

The message and the word put three factorisations of one matrix on two axes and found they rank differently on each: 48, 4 and 4 communication rounds against 1,170, 1,170 and 2,160 words. The method with the fewest rounds — Cholesky QR, which forms G = AᵀA in a single all-reduce and factorises it locally — pays for that with the site’s oldest currency.

κ(AᵀA) = κ(A)². The implied orthogonality of Cholesky QR’s Q has a fitted slope of 1.95 in κ where the Householder sweep and the reduction tree stay at 1.00, and past κ ≈ 3·10¹⁰ the Gram matrix is numerically indefinite and the factorisation stops.

The repair is one line.

The one line #

Q₁R₁ = A by Cholesky QR. Then Q₂R₂ = Q₁ by Cholesky QR again. Take Q = Q₂ and R = R₂R₁.

Fitted slopes: 1.95 for one pass, 0.93 for two. The square is gone, and what is left is κu, which is what a backward-stable factorisation costs and what the other two methods have.

The second column is not merely improved. At every condition number in the table it is at or below the Householder sweep’s — the sweep being the reference implementation of a backward-stable QR, and the method Cholesky QR was supposed to be a cheap and inferior substitute for.

Why a second pass is enough, and a third is not #

The mechanism is short enough to state and is the reason this is not a trick.

The first pass returns a Q₁ whose orthogonality error is of order κ²u. That means κ(Q₁) is about 1 + κ²u rather than κ: it is nearly orthogonal, and how nearly is exactly what the first pass got wrong. The second Cholesky QR squares that condition number — and squaring a number that is nearly one leaves a number that is nearly one.

So the second pass has an easy problem precisely because the first pass had a hard one, and the composition lands at u rather than at κ⁴u or at anything else the arithmetic might have suggested.

The two condition numbers, measured rather than asserted:

The first two columns agree to three digits over twelve decades, which is the squaring happening exactly as the algebra says — and the last column is the reason the second pass is cheap. The matrix handed to it has a condition number of one, to four decimal places, at every κ where the first pass returns anything usable at all. At κ = 10⁸, where the first pass’s Q is 0.37 from orthogonal, κ(Q₁) is 1.26.

Squaring 1.26 is 1.59. Squaring 10⁸ is 10¹⁶, which is past 1/u and is where the method stops.

A third pass has nothing left to do, and that is the refusal this essay’s library publishes: feeding the build the claim that a third pass improves the accumulated orthogonality by an order of magnitude must fail, because the second pass has already reached the floor a backward-stable factorisation has.

What it costs, counted #

One more all-reduce of n² numbers.

At p = 8 and n = 8. The sweep needs one reduction per column, so its rounds grow with n; the two Cholesky runs need one and two whatever n is.

A quarter of the sweep’s rounds for the same orthogonality is the headline, and the honest sentence beside it is that the reduction tree spends 3 rounds and 252 words and reaches the same place. On both counts this file measures, the tree is the better method.

What Cholesky QR twice has instead is its shape: two Gram matrices, two Cholesky factorisations and two triangular solves, every one of them a dense matrix–matrix operation, with no tree, no user-defined reduction operator and no dependency between processors beyond two collective calls that every library already implements. That is an argument about what a machine is fast at rather than about what it counts, and this essay does not make it — it counts.

The boundary is the first pass’s boundary #

Running it twice cannot repair a first pass that refused, and this is where the honest reporting has to be careful, because the first pass’s failure is not a threshold.

It refuses at 3·10¹⁰ and 10¹¹ and returns an answer at 3·10¹¹ and 10¹² — an answer on a worse-conditioned matrix than one it declined, which is the erratic boundary the communication field already recorded for one pass and which two passes inherit unchanged.

And where it does return past the cliff, the second pass is better than the Householder sweep: 1.2·10⁻⁶ against 3.5·10⁻⁵ at κ = 3·10¹¹. The first pass’s Q there is 1.62 from orthogonal — worse than useless — and the second pass turns that into an answer an order of magnitude better than the reference method’s.

That combination is worth sitting with. The method either declines or returns something good, and which of the two it does is not predictable from κ. A code that treats a refusal as a signal to fall back is doing the right thing; one that treats κ as a guide to whether to try is not.

Every method reconstructs, including the one whose Q is worthless #

The site’s rule earns its place here more than anywhere.

At κ = 10⁸ all four methods reconstruct A to about 10⁻¹⁶:

The residual — the quantity anybody checks a factorisation with — is identical to two significant figures across a method whose Q is 37% wrong and three whose Q is right to nine digits.

That is the refutation this essay publishes, and it is not a subtle one: a residual says nothing about orthogonality. The reason is structural rather than numerical. Q = AR⁻¹ by construction, so QR = A whatever R is, as long as R is invertible; the residual is testing the triangular solve and nothing else.

Two passes, and where the idea comes from #

The construction has a name and a literature and this essay measures it rather than deriving it, so it is worth saying what was already known and what was checked.

Known: that Cholesky QR’s orthogonality degrades like κ², that a second pass restores it while κ²u < 1, and that the composition is what the published analysis calls CholeskyQR2. None of that is a finding here.

Checked: the two slopes, on this site’s own matrices, with the site’s own implied-orthogonality measure — 1.95 and 0.93 rather than 2 and 1. The claim that the second pass reaches the Householder sweep’s accuracy rather than merely improving, which the table above makes precise and which is stronger than the usual statement. The round and word counts, from the same network model the communication field uses for the other three methods, so the four are comparable. The behaviour past the cliff, which is where the published analysis stops and where the erratic refusals live. And the residual, which is identical across a method whose Q is 37% wrong and three whose Q is right.

The last of those is the one this site adds. Every published comparison of these methods reports orthogonality; the reason the residual belongs beside it is that the residual is what a reader would otherwise check, and it says nothing at all.

What the second pass does not fix #

Three things, and the third is the one that would matter in practice.

It does not make the residual mean more. Both passes produce a Q from a triangular solve, so QR = A either way and the reconstruction is uninformative either way. A code that wants to know whether its Q is orthogonal has to measure ‖QᵀQ − I‖, which costs a product the factorisation was trying to avoid.

It does not extend the range. The cliff is where the first pass’s Gram matrix goes indefinite, and that is a property of κ(A)²u alone. Two passes have exactly the range one pass has.

And it does not help a matrix that is rank-deficient rather than ill-conditioned. Everything here is measured on matrices with a prescribed condition number and full numerical rank. A tall-skinny matrix with a genuinely dependent column has a singular Gram matrix, the first Cholesky refuses at that column, and there is nothing for a second pass to be run on — where the reduction tree and the sweep both return an R with a zero on the diagonal and let the caller decide what that means.

That last one is the practical reason the tree stays the default in libraries, and it is worth stating because it is not a numerical argument at all: it is about what a routine does when the input is not what it was promised.

Where this sits among the site’s second passes #

Running something twice to buy back accuracy is a shape this site has met three times now, and the three are worth laying beside each other because they are not the same trade.

Iterative refinement solves the system again against the residual of the first solve, and buys back the digits a low-precision factorisation lost — but only while κu < 1, and the mixed-precision field measured the threshold at exactly that boundary.

Reorthogonalising a Gram–Schmidt step projects a second time against the vectors already orthogonalised against, and takes ‖QᵀQ − I‖ from 10⁻² to 10⁻¹⁶ — the twice is enough, measured table the standard pass added to the orthogonality field.

And this. A second Cholesky QR, which is the same idea one level up: the object being re-processed is a whole factorisation rather than a vector or a right-hand side.

All three work for one reason. The first pass produces something whose error is a known function of the input’s conditioning, and the second pass sees an input whose conditioning is that error rather than the original problem’s. Composing a squaring with a nearly-one leaves a nearly-one — which is why the answer is always twice, in all three cases, and never three times.

When to use which #

The measurements support a short table rather than a recommendation, and the entries are about what a machine is slow at rather than about the matrices.

If the reduction is cheap — a small machine, a fast network, a factorisation that is not on the critical path — the Householder sweep is the reference and there is no reason to leave it. It costs n reductions and nothing else on this page beats it for simplicity.

If the reduction is expensive and κ is under about 10⁷, two Cholesky passes give the sweep’s orthogonality for two collective calls, in an arithmetic that is entirely matrix–matrix.

If κ is unknown, the refusal is the feature. The first pass declines rather than lying — usually — and a code that falls back on a refusal is correct; a code that trusts κ to predict the refusal is not, because the boundary is erratic and it returned an answer at 3·10¹¹ having refused at 10¹¹.

And if the reduction is expensive and the matrix might be badly conditioned, the reduction tree does everything this method does with fewer rounds and fewer words, and it never refuses. That is the conclusion this essay’s own counts point to, and the argument against it is one about implementation shape rather than about arithmetic — which is the honest place to leave it.

What the residual is measuring instead #

One more observation about the residual table, because it explains why every number in that column is the same.

Q is formed as AR⁻¹ in all four methods — that is what impliedOrthogonality does, and it is what a caller who solves a least-squares problem through R is implicitly doing. So QR is A R⁻¹ R, and the residual ‖A − QR‖ is measuring the accuracy of a triangular solve followed by a triangular multiply, which is a backward-stable pair whatever R happens to be.

That makes the residual a test of the triangular arithmetic and nothing else. It would catch an R with a garbage entry, a wrong dimension, a transposed factor. It cannot catch an R that is a perfectly good triangular factor of the wrong Gram matrix, which is exactly what Cholesky QR returns past its cliff.

The measurement that does catch it is ‖QᵀQ − I‖, and it costs a product the whole family of methods was designed to avoid — n² inner products of length m, which is the same order as the factorisation itself. That is the reason the check is not routinely made, and the reason this site prints it on every figure that draws a factorisation.

What the drag does #

The slider is the number of columns, and it moves one thing and leaves another alone.

What moves is the sweep’s cost: one reduction per column, so at n = 16 the Householder sweep spends 48 rounds where the two Cholesky passes still spend 6. The advantage the second pass has in synchronisations grows linearly with the width.

What does not move is either slope. κ² is a property of forming AᵀA and κ is a property of not having to, and neither depends on the shape of the matrix — which is why the figure asserts the two slopes at every position of the slider rather than at one.

The slider starts at six columns rather than four, and the reason is the cliff. At four columns the Gram matrix of a κ = 10⁸ matrix is already numerically indefinite, the first pass refuses, and the curve cannot be drawn through a refusal — which is the same boundary the table above is about, arriving in the figure’s own assertion.