A problem with no answer

A pencil has n eigenvalues because det(A − λB) is a polynomial of degree n. The previous essay was about the degree dropping when B is singular, and about the roots that go to infinity when it does.

The previous essay was about the degree dropping. There is one more thing that polynomial can do, and it is worse than dropping degree. It can be identically zero.

That happens as soon as A and B have a common null vector x: then (A − λB)x = 0 for every λ at once, so A − λB is singular for every λ, so the determinant is zero for every λ. Every complex number is an eigenvalue. The set of eigenvalues is the whole plane, which is another way of saying the question has no content.

Such a pencil is called singular — a word doing a different job from the one it does about a matrix, and one worth keeping separate: a pencil with a singular B is a perfectly good problem, and a singular pencil is not a problem at all.

Where it comes from, which is a modelling error with no symptom #

Nobody sets out to build one. Three ordinary things produce them.

A redundant equation. Two rows of a model that say the same thing — the exact-duplicate case a factorisation with nothing to pivot for meets from the other side, because a constraint was written twice or because a conservation law was imposed both globally and cell by cell. The two rows are dependent in both matrices at once, which is exactly the condition.

A variable nothing depends on. A state that appears in no equation, usually because a term was dropped during simplification or a coupling was set to zero to test something and never restored. Its column is zero in both matrices.

Or a symmetry that has not been quotiented out. A structure with a free rigid-body motion, a circuit with no ground node, a Laplacian with no boundary condition. The constant vector is in the null space of both operators, and it is there for a physical reason.

The last is the common one, and it is uncomfortable because the model is right. A floating structure really does have six rigid-body modes; the model really is invariant under them; and the pencil really is singular until somebody pins it down. What is wrong is asking an eigensolver about it.

What the exact route says, which is nothing at all #

The polynomial is computed here in BigInt rationals, by evaluating the determinant at n + 1 integer points and interpolating, with nothing rounded anywhere.

For the pencil in the hero it comes back with zero coefficients. Not small coefficients. Not coefficients below a threshold somebody chose. The list is empty, because every value the determinant was evaluated at was exactly zero and the interpolation of a function that is zero at n + 1 points is the zero polynomial.

That is a fact rather than a judgement, and it is the only unambiguous statement anybody is going to get about this pencil. Everything else in this essay is what happens when it is not available.

What a solver says, which is six numbers #

Perturb A and B by ε times a random matrix. The pencil is now regular with probability one — the common null vector is gone, det(A − λB) has degree n again — and every routine that is asked returns n eigenvalues.

At ε = 10⁻⁸ the residuals of those eigenvalues are at most 1.8·10⁻⁹.

That number is the whole problem. It is small by any standard a caller applies. It is what a backward-stable solver produces on a well-posed problem. There is nothing about it that says the computation is meaningless, and every check a careful user would run comes back clean.

Then change the seed.

The eight rows in the hero are eight seeds. No two agree. The values are spread over forty-four — which is larger than any of them is — and the spread is 2.4·10¹⁰ times the residual.

Sweeping ε does something worse than leaving the scatter unchanged: it moves the two quantities in opposite directions.

The residual tracks ε and the spread does not. Worst residual 1.39·10⁻⁷ at 10⁻⁶ against 1.8·10⁻⁹ at 10⁻⁸ — two orders for two decades — while the eight returned values stay scattered over 61.5 and then 66.4. So the ratio between them, which is what any accuracy claim here would rest on, grows as the perturbation shrinks rather than settling.

Across ε = 10⁻⁶, 10⁻⁸, 10⁻¹⁰ and 10⁻¹² the worst residual falls 1.39·10⁻⁷, 1.8·10⁻⁹, 1.76·10⁻¹¹, 1.48·10⁻¹³ — cleanly, by two decades a step, exactly as a backward-stable method should. Over the same range the spread of the answers rises: 61.5, 66.4, 6,170 and 4.6·10⁶.

So refining the model makes the report better and the answer worse, monotonically. That is the most hostile shape a diagnostic can have. A caller who suspects trouble at ε = 10⁻⁶ and responds by tightening their perturbation — the obvious move, and the right one on every well-posed problem this site has measured — will watch the residual improve by six orders while the answers they are comparing move apart by five. Every instinct that works elsewhere points the wrong way here, and the ratio of spread to residual runs 4.4·10⁸, 3.7·10¹⁰, 3.5·10¹⁴, 3.1·10¹⁹ across the four: eleven orders of magnitude of divergence between what is measured and what is true.

The exact route says the same thing at every stop and says it in one integer: the degree of det(A − λB) comes back as −1, identically zero, at ε = 10⁻⁶ and at 10⁻¹² alike. The perturbation regularises the pencil numerically and never touches the fact that it is singular.

The scatter does not converge #

That last pair is the difference between this and an ill-conditioned problem, and it is worth being precise about.

An ill-conditioned regular pencil has eigenvalues that move a great deal under a perturbation, and they move towards something: halve ε and the answers get closer to a limit, at whatever rate the conditioning allows. The answer exists and the computation is having difficulty reaching it.

Here there is no limit. Across six decades of ε the spread does not shrink. The answers are rearranged by the seed, at every size of perturbation, because there is no answer for them to be approximations to — the perturbation is not revealing a hidden true spectrum, it is creating one, and it creates a different one each time.

This is what an ill-posed problem looks like from inside a solver, and this collection has met the shape before in the regularisation field: a problem whose data does not determine its answer, where something outside the data has to choose. The difference is that a regularised problem knows it is choosing. Here nothing knows.

The three things a caller checks, and what each of them says here #

Put the ordinary diligence beside the ordinary output and the gap is the point of the essay.

Did it converge? Yes. The eigensolver reached its tolerance in the usual number of iterations, with no warning and no fallback path taken.

Is the residual small? Yes — below 1.8·10⁻⁹, which is what a backward-stable computation on a perturbed matrix produces and is several orders better than most people demand. A small residual is not a small error is the general form; this is its extreme case.

Are the eigenvalues plausible? Yes. They are ordinary real numbers spread over a range a modeller would look at without concern. Nothing is at 10³⁰⁰, nothing is NaN, nothing is complex when it should not be.

Three checks, three passes, and the answer is noise. The only check that fails is the one nobody runs: run it again with a different seed. That is why the last section of this essay is about a habit rather than about an algorithm.

Why the residual cannot see it #

The residual asks whether the computed pair (λ, x) nearly satisfies Ax = λBx on the perturbed matrices, and it does, because the perturbed pencil is a genuine regular pencil with genuine eigenvalues and the solver found them accurately.

Every word of that is true and none of it is about the pencil that was handed over.

This is the sharpest instance this collection has of a statement it has been making since its second essay: a small residual is a statement about the problem the algorithm solved, not about the problem was posed, and when the two differ by a perturbation of 10⁻⁸ that is normally an excellent bargain. Here it is not a bargain at all, because the map from the problem to the answer is not continuous, so being close to the problem buys nothing about being close to the answer.

What can be detected, and what it costs #

The exact polynomial is not available on real data, so the question is what a float computation can say. Three things, in increasing order of cost and usefulness.

The cheap check: a common null vector. Compute the smallest singular value of A and of B, and if both are near zero, check whether their null directions agree. That catches the construction here and every case where the redundancy is a single vector. It costs two singular value decompositions and it misses everything more subtle than one shared direction.

The honest check: the smallest singular value of the pencil over a set of λ. A pencil is singular exactly when σ_min(A − λB) is zero for every λ. Sampling that function over a grid and finding it uniformly tiny is strong evidence; finding it small at isolated points is an ordinary spectrum, which is a decision about a zero taken over a region rather than at a point. It is the same idea as a pseudospectrum, which this collection has already drawn, and it costs one singular value decomposition per grid point.

And the right check, which nobody runs: the staircase. There is a canonical form for pencils — the Kronecker canonical form — that separates the regular part from the singular part exactly, and the staircase algorithms that compute it are the correct answer to this problem. They are also delicate, they need a sequence of rank decisions each of which can be wrong, and they are not what a caller gets by asking for eigenvalues.

The practical position is therefore uncomfortable and worth stating plainly: the standard tool returns an answer, the answer is meaningless, the diagnosis is available, and it is not run by default.

The perturbation is doing the choosing #

There is a sharper way to say what the perturbation does, and it makes the failure sound less like bad luck and more like what it is.

A singular pencil is a point on a surface in the space of all pencils, and the surface has codimension greater than one — moving off it in any direction produces a regular pencil, and the eigenvalues that come out depend on the direction rather than on the distance. Halving ε moves half as far off the surface along the same ray and returns the same spectrum; changing the seed changes the ray and returns a different one.

So the answers are a function of the perturbation’s direction, which is exactly the thing a rounding error has no opinion about. Every arithmetic path through the computation is a different direction, which is why the same solver on the same data can give different answers on two machines, in two compilers, or with two orderings of the same sum — and the differences are the size of the answers rather than the size of the rounding.

A reduction that changes the order of a sum changes the result, and that essay’s differences are in the last few digits. Here the same mechanism produces differences in the first — which is why two machines disagree writ large enough to change a conclusion rather than a digit.

Two seeds are the cheapest thing available #

If nothing else is done, one habit costs almost nothing and catches this.

Run it twice. Perturb the pencil by a different tiny random amount and compare the spectra. On a regular pencil, however ill-conditioned, the two runs agree to whatever the conditioning allows and the disagreement is bounded by it. On a singular pencil they disagree by the size of the answers.

That is this site’s two-routes habit applied to a case where the second route is the same route with a different seed — which is the weakest form of it, and here it is enough, because the failure being detected is precisely a failure of the answer to be a function of the input.

It also has the property that makes a check worth having: it produces a number rather than a verdict. The spread across seeds, divided by the size of the eigenvalues, is a continuous measure — and what that number actually is turns out to be the least satisfactory part of the whole essay, which the next section measures.

The three detectors, ranked by reliability rather than by cost #

The list above is ordered by cost, and the implication is that cost buys reliability: the cheap common-null-vector check “misses everything more subtle than one shared direction”, the grid check is “strong evidence”, the staircase is “the right check”. Running the first, the second and the seed spread on the hero’s own pencil across fourteen decades of perturbation puts them in a different order.

The cheap check is the one that never fails. The two null directions agree to seven decimal places out to a perturbation of 10⁻⁴ and to three at 10⁻², where the perturbation is a per cent of the matrix and any modeller would call the data ruined. Two singular value decompositions, no grid, no canonical form, and it holds across the whole range.

The grid check tracks it and costs sixty-one decompositions rather than two, so on this pencil it buys nothing the cheap one did not already have. Its value is the one the essay states — it catches redundancy that is not a single shared direction — and this construction has a single shared direction, so the measurement is fair to the check rather than to the claim.

And the seed spread, the one this essay goes on to recommend, is the noisiest of the three. It does fire: every value is thirty or more, against something near the rounding for a regular pencil. But it ranges from 34 to 7.1·10⁹ across the sweep, so “around one here” is not a calibration — it is wrong by between one and ten orders depending on a perturbation size the caller chose for unrelated reasons. The check answers is this pencil singular correctly and answers how singular with a number that means nothing.

That reorders the practical advice. Compute σ_min of A and of B and compare their null directions, always — it is two decompositions on a problem that is about to get an eigensolve, so it is free by comparison, and on the redundancy that actually arises it is decisive across fourteen decades. Run the two seeds as well, because it costs one extra solve and catches redundancies of more than one dimension that the cheap check will miss. And read its output as a yes or a no rather than as a measure of anything, because the number it produces is a function of the perturbation and not of the pencil.

What a library actually does about this #

It is worth being fair to the software, because the situation is better than “nobody checks” suggests and worse than the check being present would imply.

LAPACK’s generalised eigenvalue routines return, alongside each eigenvalue, a pair of numbers that bound its conditioning — a reciprocal condition number for the eigenvalue and one for the associated eigenvector, available on request. On a singular pencil those come back tiny, because the eigenvalue is infinitely sensitive and the estimate says so as loudly as a finite number can.

So the diagnosis is available from the standard tool, at the cost of asking for it.

What it is not is the default. The routine most callers reach for returns the eigenvalues and nothing else; the conditioning estimates are behind a different entry point with more arguments; and the documentation that explains what a tiny one means is not the documentation somebody reads while getting a model to run.

That is a recognisable shape and this collection has met it twice already. A condition estimate that is computed and not printed; a residual that is available and not asked for; a backward error that every library could report and none does by default. The information exists, costs almost nothing, and sits one function call away from the person who needs it.

The refusal #

The assertion is fed a request for the eigenvalues of a singular pencil.

That is the failure this essay is about, in the one place where it can be caught: the routine has the exact characteristic polynomial in front of it, has found that it has no coefficients, and is being asked to return a finite set of roots anyway. Refusing is the difference between a routine that answers a question and one that answers a different question quietly.

The refusal is worth having even though the exact polynomial is rarely available, because it settles what the routine’s contract is. A function that returns eigenvalues of a pencil is claiming the pencil has some, and a function that cannot make that claim should say so rather than return the roots of whatever it was handed.

What is next #

Three essays about a pencil, and the thread running through them has been what a zero in the second matrix means: a constraint, an eigenvalue at infinity, or the whole question falling over. The next pair changes subject and keeps the theme. They are about a claim this site has been making in a source comment since its first month — that one particular algorithm computes small singular values to high relative accuracy — which turns out to be true, and to be true of a class of matrices rather than of the algorithm.