A test with no answer in it

Almost every measurement in this collection is a comparison against something known: a matrix built with a prescribed spectrum so that its condition number is an axis rather than an accident, an inverse available in closed form so that the true answer can be written down rather than approximated, a second computation carried at a higher precision and treated as truth. Each is a way of manufacturing a right answer so that a computed one can be subtracted from it, and none of it is available to the person calling a library. Somebody holding a routine and a matrix that arrived from a measurement has no reference answer at all, and the usual conclusion is that nothing can be checked except the residual.

There is one thing such a caller can still do, and it needs no answer and no second precision: ask the routine the same question twice, in two forms the problem cannot tell apart, and check that the replies agree. The question here is orthogonalise this, and the two forms are a matrix and the same matrix with its columns written in the opposite order. Nothing about a nearest orthogonal matrix depends on which column was written first, so a routine whose answer does depend on it is answering a different question — and no reference answer was needed to establish that.

Run that way the check does not merely rank the two standard orthogonalisations. It separates them completely, at every conditioning drawn.

The identity is one line, and the test is one extra call #

The upper grid is zero not because the polar factor is computed carefully but because it is an identity, and the identity falls out of the singular value decomposition in a line.

Write A = PΣVᵀ. Permuting the columns of A by an orthogonal permutation matrix Π — the matrix the figures label P, since they have no singular vectors to name — gives AΠ = PΣ(VᵀΠ) = PΣ(ΠᵀV)ᵀ, which is a singular value decomposition of AΠ with the same P, the same Σ, and V replaced by ΠᵀV. The orthogonal polar factor is the product of the first and third of those, so it comes back as P(ΠᵀV)ᵀ = PVᵀΠ, which is the original factor with the same permutation applied. The singular values never moved; only the right basis was relabelled, and the factor inherits the relabelling exactly.

Nothing of the sort is available to a QR. Householder introduces zeros one column at a time and Gram–Schmidt projects each column against the ones already finished, so both leave the first column of A essentially alone and subject the last to every projection in the sequence. Permuting the columns changes which column is which, and there is no algebra to recover: the two Q’s answer two different questions about two different nested sequences of subspaces. The essay that established which orthogonal matrix is nearest put that difference at nine to seventeen per cent of distance and observed in passing that the ordering dependence is the much larger fact.

What the identity buys is a test with the shape a caller needs. It costs one more run of the routine under test, a reordering of the columns and one Frobenius norm — a data movement and a sum of squares, so nothing in the overhead grows faster than n². There is no truth to construct and no second implementation to find. The alternative is to build a reference answer, which for this question means computing a singular value decomposition, which is very often the routine under test.

At κ = 10 and n = 6 the two sides read 1.53·10⁻¹⁵ and 2.353, fifteen orders of magnitude apart. No threshold has to be argued over, because there is nothing in between.

What the test reads when there is nothing to orthogonalise #

The figure above is the control, and it is where the test is weakest. A matrix with κ = 1.05 is already close to orthogonal, so every method returning something orthogonal returns something close to A and therefore close to what every other method returned. The two QR answers are then close to each other for a reason that has nothing to do with the algorithm being order-independent.

The collapse is orderly and it can be followed down. At n = 6, the reversal test reads 0.0696 at κ = 1.05, 0.0142 at κ = 1.01, 0.0014 at κ = 1.001 and 0.0001 at κ = 1.0001. As a fraction of the 2.449 those matrices measure, that is 2.8 per cent, 0.6 per cent, 0.06 per cent and 0.006 per cent. A caller who wrote the obvious check — reverse the columns, take the norm of the difference, complain if it exceeds a hundredth — passes a Householder QR on every one of those matrices.

That is not a hypothetical range. A rotation matrix that has drifted over a few thousand integration steps, a basis that has lost its orthogonality inside a subspace iteration, a transformation fitted to slightly noisy point pairs: all of them sit within a thousandth of orthogonal, and all of them are exactly the inputs on which the reversal test returns a number small enough to look like rounding. The diagnostic is quietest on the problems it was reached for.

Two things usually run together are worth separating. The margin between the identity’s residual and the violation does not collapse: at κ = 1.05 the two sides read 1.80·10⁻¹⁵ and 0.0696, thirteen orders, and at κ = 1.0001 they read 2.07·10⁻¹⁵ and 0.0001, eleven. What collapses is the absolute size of the violation, so only a caller comparing it against a fixed number is misled — and the repair is available to that reading and not to the other.

The floor is not a matter of size #

The obvious suspicion about a diagnostic that goes quiet is that it is running out of room, and at n = 3 there are only three columns to reorder. The sizes say otherwise.

‖Q‖_F is exactly √n for any orthogonal Q, which gives the readings a scale that does not have to be estimated: 1.732, 2.000, 2.236, 2.449, 2.646 and 2.828 at n from three to eight. Measured against it, the reversal test at κ = 10 reads 0.687, 0.932, 1.043, 0.961, 1.023 and 0.959 of the matrix’s own size, and at κ = 1.05 the same six sizes read 0.016, 0.025, 0.033, 0.028, 0.029 and 0.027 — a factor of about thirty-five between the two regimes at every size except the smallest, where it is forty-three.

So the collapse is a property of the conditioning and not of the dimension. Three columns are enough to expose the ordering dependence completely — 1.190 against 1.732 is more than two thirds of the matrix’s own size — and eight columns do not rescue the test at κ = 1.05, where it reads 0.077 against 2.828. Whatever the floor is, it is not a shortage of columns to permute.

The polar side reads 2.50·10⁻¹⁶ at n = 3 and 2.55·10⁻¹⁵ at n = 8 at κ = 10 — an order of magnitude for a factor of eight in the entries being summed, which is what a Frobenius norm over n² rounding errors should do.

The floor is in the units, not in the test #

The reading collapses because it is being divided by the wrong thing. ‖Q‖_F is the size of the answer, and the answer has unit-length columns whatever the input: the same 2.449 at κ = 1.0001 as at κ = 10⁴. What actually shrinks with the conditioning is the correction — the distance ‖A − Q‖ from the input to whatever the routine returned — and that is what the ambiguity should be compared against, because it is the size of the thing being computed. It costs nothing extra: a caller running the reversal test already holds A and Q, so ‖A − Q‖ is one subtraction and one sum of squares, and it needs no reference answer either.

In those units at n = 6 the reversal test reads 0.867 at κ = 1.0001, 0.867 at 1.001, 0.870 at 1.01, 0.883 at 1.05, 0.997 at 1.5, 1.319 at 10, 1.398 at 100, 1.412 at 1,000 and 1.409 at 10⁴. The floor is gone. The number a caller would threshold on varies by less than a factor of two across a range where the reading in the other units varied by a factor of thirty-five thousand.

At n = 8 the same column reads 0.846, 0.846, 0.849, 0.862, 0.973, 1.328, 1.462, 1.495 and 1.496; at n = 3 it reads 0.491 through 1.136. The smallest size sits lower and is still flat, which is the shape that matters: a diagnostic whose reading is roughly constant across the regime of interest can be given a threshold, and one whose reading spans four decades cannot.

What the constant says is the sharpest thing the test measures. The two QR answers differ from each other by very nearly as much as either differs from A, so the correction applied to a drifted rotation is matched in size by a comparable correction in an unrelated direction: the arbitrary part of the answer is the size of the answer. That reading is available to a caller with no truth at all, and it is a stronger statement than the nine to seventeen per cent of extra distance that needed a reference answer to establish.

The price of the better units is at the other end, and it is real. Dividing the polar side by ‖A − U‖ divides a residual near 2·10⁻¹⁵ by a correction that goes to zero, so the identity’s own reading rises as the matrix approaches orthogonal: 1.4·10⁻¹¹ at κ = 1.0001 against 2.5·10⁻¹⁴ at κ = 1.05 and 9.5·10⁻¹⁶ at κ = 10. The margin between the two sides is therefore ten orders at the extreme left of the range rather than thirteen. Ten orders is not a difficulty; it is the honest cost of the units that removed the floor.

The identity’s own side is not zero either #

Everything above treats the polar reading as though it were a zero. It is not, it grows, and the rate at which it grows is what decides how far the test can be pushed.

Across the slider at n = 6 the polar reading is 1.80·10⁻¹⁵, 1.53·10⁻¹⁵, 2.78·10⁻¹⁵, 5.28·10⁻¹⁵ and 7.68·10⁻¹⁴ at κ = 1.05, 10, 100, 1,000 and 10⁴: flat through three decades and then climbing. Continuing past the right-hand stop shows where it goes — 1.19·10⁻¹³ at κ = 10⁵, 1.69·10⁻¹² at 10⁶, 1.44·10⁻¹⁰ at 10⁸, 6.79·10⁻⁹ at 10¹⁰, 2.02·10⁻⁷ at 10¹² and 1.01·10⁻⁵ at 10¹⁴. Ten decades of conditioning cost about eight decades of the identity, a little below linear in κ and consistent with the decomposition delivering the small singular directions to an accuracy that degrades with the spread.

The violation, meanwhile, saturates: 2.353, 3.048, 3.325, 3.435, 3.459, 3.464, and then it stops, reading 3.4641 at κ = 10⁸ and 3.4641 at 10¹⁴. That ceiling is not arbitrary. Two independently drawn orthogonal matrices of size n satisfy ‖Q₁ − Q₂‖²_F = 2n − 2 tr(Q₁ᵀQ₂), whose expected value is 2n, so √(2n) = √12 = 3.4641 is exactly the distance at which two answers carry no relationship to each other at all. QR’s two answers reach that value and stay there, which is the strongest form the finding takes: past a certain conditioning, reversing the columns produces an answer as unrelated to the first as one drawn at random would be.

So the test’s headroom is a measurable quantity that narrows: 3.9·10¹³ at κ = 1.05, 1.5·10¹⁵ at κ = 10, 4.5·10¹³ at κ = 10⁴ and 3.4·10⁵ at κ = 10¹⁴. Five orders is still an enormous margin, and it was thirteen. The headroom is not infinite in principle and it does not run out inside the range where the input is still a matrix rather than a numerically singular one — a bound worth stating rather than assuming, since the same cannot be said of every diagnostic that reports a small number.

Which permutation is not a free choice #

The test as described says reverse the columns, and a caller might reasonably swap two of them instead. The choice is not free, and the size of the difference is the second floor in this essay.

All 720 permutations of six columns can be run at κ = 10, and the distribution is wide. Excluding the identity, which returns exactly zero as it must, the QR reading has a minimum of 0.0721, a tenth percentile of 1.1969, a median of 1.8211 and a maximum of 2.3945. The reversal reads 2.3529, within two per cent of the largest reading available — not luck, since reversal moves every column furthest from its original position and the ordering dependence of a sequential algorithm is largest exactly there.

The weak end of the distribution is the warning. Ninety-three per cent of the 719 non-trivial permutations produce a reading above 1, and the remaining seven per cent produce readings a threshold set for the reversal would pass. Swapping single adjacent columns of the same matrix returns 1.3405 for the first pair, 0.0721 for the second, 0.0738 for the third, 0.4868 for the fourth and 0.6721 for the fifth. A caller who swaps columns two and three reads 0.072, the same number the reversal returns at κ = 1.05. Two ways of running the test out of its useful range give the same misleading reading, and only one of them is about the matrix.

The polar side is untroubled by any of this. Over all 720 permutations its worst reading is 2.68·10⁻¹⁵, against 1.53·10⁻¹⁵ for the reversal alone — the identity holds for every reordering, as the derivation says it must, and the spread over the whole symmetric group is less than a factor of two.

What the test cannot see, which is most things #

An invariance test is a test of equivariance, and equivariance is not correctness. How little the check establishes on its own is worth stating bluntly, because the numbers are extreme.

A routine that returns its input unchanged passes perfectly. Its two answers are the input and the permuted input, so the difference is exactly zero — not 10⁻¹⁵, but zero, since permuting columns moves floating-point numbers without altering them. That same matrix measures 1.9694 for the quantity orthogonality actually is, ‖XᵀX − I‖, and the reversal test has nothing to say about it.

Worse, a routine that is genuinely trying passes too. Newton’s iteration for the polar factor is built from matrix products and inverses, all of which commute with a column permutation exactly, so it is equivariant at every step whether or not it has converged. Stopped after one step it returns a matrix 6.2277 from A with an orthogonality error of 26.5187; after two steps, 3.2641 from A with an error of 6.3000. Both read exactly zero on the reversal test.

The repair is the point: the reversal test is one half of a pair, and the other half is the orthogonality residual. ‖XᵀX − I‖ needs no reference answer either, costs one matrix product, and is the check every library already runs. It reads 1.43·10⁻¹⁵ for the sign-fixed QR and 2.32·10⁻¹⁵ for the polar factor and separates neither, which is precisely why the second check was needed. Each passes exactly the candidates the other rejects.

The contrast is in what each one needs. Ranking four candidates by distance requires knowing what nearest means, which requires the controls beneath the bars, which requires the ability to draw random orthogonal matrices and the willingness to spend two hundred factorisations establishing an ordering the algebra already guarantees. The invariance test asks a smaller question, needs none of that, and answers it with a margin of thirteen orders rather than the fifteen per cent separating the top two bars.

Both checks together still leave exactly one matrix #

The two checks being complementary invites the conclusion that together they pin the answer. They very nearly do, and what remains is a specific matrix at a measurable distance, which is a better thing to know than a reassurance.

Suppose a routine returns U(A)·S for some fixed orthogonal S. It is orthogonal, so the residual check passes, and its reversal reading is ‖U(A)ΠS − U(A)SΠ‖ for the reversal Π, which vanishes exactly when S commutes with Π. Every orthogonal matrix commuting with the chosen permutation is therefore an undetectable impostor, and the reversal — being its own inverse — commutes with a great many. Taking S to be the reversal itself gives a routine that reads 1.53·10⁻¹⁵ on the test and sits 2.8768 from A against the true factor’s 1.6054: 79 per cent worse, and invisible.

Running every permutation instead of one removes almost all of them, because the matrices commuting with the whole symmetric group are only those of the form αI + βJ with J the all-ones matrix. Of those, exactly two are orthogonal: the identity, and S = I − (2/n)J, the reflection in the all-ones direction. That second one is a genuine survivor, and it was measured rather than argued: it is orthogonal to 2.49·10⁻¹⁵, its reading is at most 2.68·10⁻¹⁵ over all 720 permutations of six columns, and it is a completely different matrix — ‖US − U‖_F = 2.000000000000, exactly two, independently of n, because the perturbation is (2/n)J and its Frobenius norm is (2/n)·n.

The cost of the blind spot is small and it is not zero. That impostor sits 1.8535 from A against the true factor’s 1.6054, 15.5 per cent further — about the penalty a sign-fixed QR carries, arriving from a completely different direction. The complete answer-free suite, run over every permutation of the columns, therefore certifies the answer to within one specific alternative that is one sixth worse.

Closing it costs nothing, and the fix is the generalisation the identity always permitted: polar(AV) = polar(A)V holds for any orthogonal V, not only for a permutation. With a randomly drawn orthogonal right factor the polar factor still reads between 1.17·10⁻¹⁵ and 2.97·10⁻¹⁵ over two hundred draws, while the impostor reads between 1.3201 and 2.8284 with a median of 2.7104. It is caught on the first draw, every time.

And the generalisation does not help the floor. At κ = 1.05 the QR reading over two hundred random orthogonal right factors runs from 0.0274 to 0.0659 with a median of 0.0501 — every one below the reversal’s 0.0696. A nearly orthogonal matrix has almost nothing to orthogonalise, and no choice of perturbation manufactures a signal that is not in the problem. The floor belongs to the matrix; the blind spot belonged to the choice of test, and only the second is a design decision.

What follows for anything that checks a routine #

A caller with no truth is not without a test. An identity the answer satisfies is a test, it needs no reference and no second precision, and here it costs one extra call plus work growing like n².

Report the ambiguity in units of the correction, not of the answer. ‖Q‖_F is √n whatever happened and cannot shrink; ‖A − Q‖ is the size of what the routine did. The same measurement divided the two ways spans a factor of thirty-five thousand or a factor of two over the same range, and only one of those can be given a threshold. This is the reading a scaling-invariant condition number forces on a different quantity: what a number is divided by is part of what it measures.

Pair an equivariance check with a residual check, because neither is worth much alone. A routine that does nothing passes the first exactly and fails the second by 1.97; a QR passes the second at 10⁻¹⁵ and fails the first by more than its own size. That is the instinct behind the certificates a solver could return instead of an answer: several cheap properties, each of which the wrong answer fails differently.

State the range in which a diagnostic fires. This one is quietest on exactly the inputs it was reached for, and a measurement established that rather than an argument, which is also what a tolerance chosen without measuring the machine’s own variation turns out to need.

And prefer a structural invariance to an empirical one. The identity here comes from the decomposition, holds exactly over the reals and degrades only through the arithmetic, which is why its reading is 10⁻¹⁵ and not 10⁻⁶ — a different object from a run-to-run comparison of one program with itself, where the variation has no algebraic value to be compared against and the threshold has to be built from observations.

The wider lesson is one this collection keeps arriving at from different directions. A quantity that moves under something the problem does not move under is measuring the description rather than the problem, which is the argument the determinant fails and the reason a null-space basis chosen by a one-line rule can be eight orders worse than another that is equally correct. What is new here is the direction of use: the principle normally rejects a quantity, and applied to a routine it becomes an instrument. The invariance the right answer satisfies is a certificate the wrong answer cannot forge, and it can be checked by anyone holding the matrix.

Two closing qualifications, both measured. The instrument does not distinguish Householder from Gram–Schmidt: both are sequential in the columns and both fail this test by more than their own size, while the orthogonality residual that separates them so decisively says nothing here. And it does not distinguish a converged polar iteration from an unconverged one, since the inverse-free iteration is equivariant at every step. The test answers one question, exactly, for one extra factorisation, and the discipline is in knowing which question that is.