The plane survives what its vectors do not

The gap decides the eigenvector ends with a computation that cannot be rescued. Two eigenvalues a billionth apart, a perturbation of 10⁻⁶, and the returned eigenvectors are half a radian from the true ones and are decided by the perturbation rather than by the matrix. No algorithm does better, because there is no better answer: at an exactly repeated eigenvalue every vector in the plane is an eigenvector, and near-repetition inherits the indeterminacy.

The way out is not a better algorithm. It is a better question.

On the same matrices, at the same perturbation, at the gap where the individual vectors have become meaningless, the plane they span moves by 7.6·10⁻⁸ radians — and does not move any further as the cluster tightens, because the quantity governing it is the gap to the rest of the spectrum, which was never small and never changed.

The ratio between the two, at the tightest gap drawn, is six million.

What an invariant subspace is, and why it is the right object #

A subspace S is invariant under A when AS ⊆ S: applying the matrix to any vector in S produces another vector in S. The span of any set of eigenvectors is invariant, and — this is the part that matters — the span of the eigenvectors of a cluster is invariant whether or not the individual eigenvectors within it are determined.

For a symmetric matrix with eigenvalues λ₁ ≤ … ≤ λₙ, the span of the first m eigenvectors is invariant, and it is unique as long as λₘ ≠ λₘ₊₁. That last condition is the whole distinction: the plane of the first two is unique provided the second eigenvalue differs from the third, and it says nothing about whether the first differs from the second.

So the object is well posed exactly when the gap out of the cluster is nonzero, and the object that is ill posed — a particular vector inside it — is ill posed exactly when a gap within the cluster is small. Two different gaps, two different questions, and the standard way of asking (give me the eigenvectors) is the one that couples the answer to the wrong one.

Measuring how far a subspace moved #

An angle between two vectors is an arccosine. An angle between two subspaces needs a definition, and the standard one is principal angles: for orthonormal bases Q₁ and Q₂ of the two subspaces, the principal angles are the arccosines of the singular values of Q₁ᵀQ₂, smallest first, and the usual single number is the largest of them.

The definition is basis-independent, which is the property it exists for and the one worth checking before believing anything computed with it. Three checks run: a plane against itself gives zero; a plane against a rotated basis for the same plane gives zero to 10⁻¹²; and two planes sharing exactly one direction give one principal angle of zero and one of exactly π/2. Without those, a basis-dependent “angle” would produce plausible plots of nothing.

The definition also cannot be computed the way it is written, which the previous essay measures: the arccosine of a singular value near one loses half its digits, so the obvious route floors at √(2u) = 1.49·10⁻⁸ regardless of the angle. Every plane angle in this essay is below that floor. They are computed by the sine formulation — the singular values of (I − Q₁Q₁ᵀ)Q₂ — and the figure of the floor is in the previous essay because the instrument had to be established before the measurement it takes.

The measurement #

Four gaps, one perturbation of 10⁻⁶, the same runs the eigenvector table came from:

The plane’s column is constant to three digits across four decades of gap. Its bound is ‖E‖ divided by the gap to the rest of the spectrum — which is 2.99 here, the distance from 1 + gap to the next eigenvalue at 4 — and that denominator does not change when the cluster tightens, because the cluster tightening is a statement about the inside of the plane.

At the bottom row the two quantities differ by a factor of 5.7 million and the eigenvalue’s own shift, 3.87·10⁻⁷, sits between them.

That is the practical content of the entire subject. Ask for an invariant subspace rather than for eigenvectors, and a hopeless computation becomes a well-conditioned one — with no change to the matrix, the perturbation, the precision or the algorithm. What changed is which function of the output was asked for.

Where the outer gap comes into it #

The derivation is the one the previous essay sets up, carried one step further, and it is worth carrying because the answer falls straight out.

To first order the perturbation of eigenvector xᵢ is a sum over the other eigenvectors,

δxᵢ  =  Σ over j ≠ i of  xⱼ (xⱼᵀ E xᵢ) / (λᵢ − λⱼ)

with each term divided by its own gap. For an eigenvector inside a tight cluster the term whose j is the cluster partner has a denominator of 10⁻⁹ and dominates everything else by six orders of magnitude, which is why the vector swings.

Now ask what happened to the plane rather than to the vector. The plane’s movement is what is left after the components of δx₁ and δx₂ that lie inside span(x₁, x₂) are discarded — because a displacement within the plane rotates the basis and does not move the plane. And the dominating term is exactly one of those: δx₁’s largest component is along x₂, which is in the plane.

So the sum that governs the plane is the same sum with the intra-cluster terms deleted:

Σ over j outside the cluster of  xⱼ (xⱼᵀ E xᵢ) / (λᵢ − λⱼ)

every denominator of which is a gap to something outside — 2.99 at the smallest, here. The tight gap has been removed from the expression rather than bounded, which is why the plane’s angle is not merely smaller but independent of it, flat to three digits across four decades.

That is also the derivation’s version of the sentence about symmetry. The ambiguity at an exactly repeated eigenvalue is a rotation within the plane; the quantity that ignores rotations within the plane is the quantity whose expansion has no small denominator in it. The two statements are the same statement.

The object underneath: a projector #

There is a way of saying all this that removes the word “basis” entirely, and it is the form the theory is usually stated in.

The spectral projector onto a cluster is P = Σ xᵢxᵢᵀ over the cluster’s eigenvectors. It is a symmetric matrix, it satisfies P² = P, and — the point — it does not depend on which orthonormal basis of the plane the xᵢ happen to be. Rotate the two eigenvectors within their plane and every individual entry of every xᵢ changes; the sum x₁x₁ᵀ + x₂x₂ᵀ does not change at all.

So the projector is the invariant object, and the plane is just a way of describing it. The largest principal angle between two subspaces is related to their projectors by

‖P̂ − P‖₂  =  sin θmax

so the angle this essay measures is a distance between two matrices, in a norm, with no basis anywhere in it. That is the form in which “the subspace is well conditioned” is a statement about computable objects rather than about geometry.

It also explains a practice that otherwise looks like carelessness. A routine that computes a cluster’s eigenvectors and then orthogonalises them is not improving the vectors — they are as indeterminate afterwards as before. It is producing an orthonormal basis, from which the projector can be formed correctly, and the projector is what everything downstream of it actually uses: a spectral filter, a deflation, a low-rank correction, a block of a similarity transformation. None of those needs to know which basis was chosen.

Why every serious eigensolver already does this #

The result is not obscure and its consequences are visible in software, where they usually read as implementation details.

An eigensolver that reports clustered eigenvalues together, or that returns an orthonormal basis for a cluster rather than individual vectors, or whose documentation says the vectors within a cluster are “not individually determined”, is applying this. So is a solver that computes the eigenvectors of a cluster by inverse iteration and then explicitly orthogonalises them against each other: the orthogonalisation does not make the individual vectors more accurate, and it does make the basis an orthonormal basis for the right plane, which is the accurate object.

The pattern generalises past eigenvectors. Rank is a decision makes the same move for singular values: which singular vectors correspond to which singular value is not determined when two singular values are close, and the subspace spanned by the leading k is, which is why a low-rank approximation is well conditioned where its factors are not. And the best approximation there is states the Eckart–Young result about the approximation rather than about the factors, for the same reason.

The rule underneath all three: a quantity that is invariant to the ambiguity is the quantity to compute. Where a problem has a symmetry — any rotation within the plane, any sign flip, any reordering within a tie — the well-conditioned outputs are the ones that symmetry does not move.

Where the subspace stops being safe #

The plane is governed by the gap out of the cluster, so the plane is in trouble when that gap is small — which is not a special case, it is the same phenomenon one level up.

Take a spectrum of 1, 1 + δ, 1 + δ + η, and the rest far away. The plane of the first two is governed by η. The three-dimensional subspace of all three is governed by the distance to the fourth eigenvalue, which is large. So the well-conditioned object depends on where the clusters are, and a routine that reports “the leading two-dimensional invariant subspace” is making a claim that needs the gap after the second eigenvalue to be checked.

There is no way round this and no need for one: the gaps are available for nothing once the eigenvalues have been computed, and choosing which subspace to report is a decision a caller can make with them. What is not available is a residual test that would catch a bad choice, for the reason the previous essay gives — every vector in a cluster’s plane has a small residual, and so does every subspace containing it.

What to ask for, in practice #

The essay’s advice is short and it is worth separating from the theory, because a reader who takes away only “subspaces are better conditioned” will apply it in the wrong places.

Compute the eigenvalues and look at the gaps first. They cost nothing once the decomposition is done, and they are the only thing that says which of the outputs about to be used is well posed. A cluster is a set of eigenvalues whose internal gaps are small compared with u‖A‖ — or, more usefully, small compared with whatever uncertainty the matrix’s entries have, which in an application is usually far larger than the rounding.

Group the clusters and use bases, not vectors. For each cluster, keep the orthonormal basis the solver returned and use it as a block. Anything computed from the block — a projector, a filtered matrix function, a deflation, a subspace to restart an iteration from — is well conditioned. Anything that names an individual column of it is not.

Do not orthogonalise and then trust the columns. Orthogonalising a cluster’s vectors is the right thing to do and it improves the basis, not the vectors. A caller who reads column 3 of the orthogonalised block and treats it as “the third eigenvector” has undone the whole benefit, and nothing in the return value marks the difference.

And do not expect a residual to warn. Every vector in a cluster’s plane has a residual at rounding, so the convergence test the solver used cannot distinguish them and will not. That is the one place where this subject is genuinely a trap rather than a subtlety: the usual defence is available, is passed, and is measuring the eigenvalue.

The uncomfortable corollary is that a great deal of code does name individual columns, because the interface it is calling offers them and because the documentation’s remark about clusters is easy to read past. The measurements in these two essays are what that costs at a gap of 10⁻⁸: an answer half a radian from the one the caller thinks it has, with no diagnostic anywhere in the computation.

What is asserted, and what refuses #

Principal angles are basis-independent, checked three ways before anything else here is believed.

Every plane angle is inside its bound — ‖E‖ divided by the gap to the rest of the spectrum — at every gap on the sweep.

The plane’s angle does not move with the cluster’s gap, to within a factor of three across four decades, asserted alongside the eigenvector angle’s growth so the contrast is a single check.

The ratio exceeds 10⁵ at the tightest gap, asserted as a number rather than left as a description.

At an exactly repeated eigenvalue the plane is identical across three perturbations, to 1.8·10⁻¹⁵, while the individual vectors differ by up to a radian — and every one of those vectors satisfies Ax = λx to 5·10⁻¹⁵.

And the eigenvector’s gap used to bound the subspace is refused. The claim “the plane moves by as much as its vectors do” is fed to rejects() at a gap of 10⁻⁸ and fails, which is the assertion that keeps the two bounds from being quietly interchanged — the mistake that makes a well-posed computation look hopeless.

What is left #

The unsymmetric case. Invariant subspaces of an unsymmetric matrix are governed by the separation of two operators rather than by a difference of eigenvalues, which is a genuinely different quantity and reduces to the gap only in the symmetric case. Everything here is symmetric, deliberately.

The algorithms. Nothing in this essay computes a subspace by a method designed for the purpose — subspace iteration, block Lanczos, the ones that return a basis directly and never form an individual eigenvector. What is measured is the conditioning, which is the ceiling those algorithms run into.

How to choose the cluster. Deciding which eigenvalues to group is a decision about a gap, and it is the same decision, with the same character, as choosing a rank — which this site has already argued is a judgement rather than a computation. That the same judgement appears here, about the same kind of quantity, is a connection worth an essay of its own rather than a paragraph at the end of this one.