A definite pivot sits at the mean

The push a pivot needs belongs to the matrix asked whether an indefinite factorisation’s most negative pivot could serve as a shift for inverse iteration towards the most negative eigenvalue. On thirty random symmetric matrices it could not, under either of two pivoting rules: the push each matrix needed ran from 0.59 to 2.27 and no fixed push served. What did work needed no knowledge of the push at all — multiply the shift by 1.25 until the factorisation of $A - \sigma I$ has no negative pivot, which Sylvester’s law of inertia says happens exactly when $\sigma$ is below the spectrum. It landed every time, in two or three factorisations.

Its last section turned the question round. On a positive-definite matrix the factorisation is Cholesky’s, every pivot is positive, and the smallest pivot is the natural candidate for a shift towards the smallest eigenvalue. Half the answer was given in advance by an argument: each pivot is one over a diagonal entry of the inverse of some trailing block, every such entry is at most one over that block’s smallest eigenvalue, and interlacing puts that at or above the whole matrix’s smallest eigenvalue. So the smallest pivot is never below the spectrum, and the search has to push it down at least once. The other half was a prediction with a sign: “on thirty random definite matrices of the same order it needs a median of four pushes, more than either rule needs here, because a definite matrix’s pivots sit near its diagonal while its smallest eigenvalue sits far below it — and the pivot is a usable shift unpushed on none of them.”

The argument holds on every matrix. The count is wrong by a factor of two, and the reason given for it is close but not quite the reason.

Three definite families #

The earlier measurement’s matrices were random symmetric matrices of order 24 with Gaussian entries, which are indefinite. Three definite families of the same order stand in for them here, thirty matrices each:

• the same random symmetric matrices, shifted so their smallest eigenvalue is exactly one half — the nearest definite relative of the earlier control;
• Gram matrices $B^{\mathsf T}B/m$ of Gaussian $B$ with $m = 48$ rows, the sample covariance of forty-eight observations of twenty-four variables, whose spectrum fills an interval from about a twelfth to about three;
• the same with $m = 30$ rows, where the bottom of the spectrum comes much closer to zero.

Each is factorised four ways. When symmetry is not enough is why there are four rather than one: on an indefinite matrix a factorisation restricted to the diagonal can fail outright, and the symmetric rules exist to avoid that. On a definite matrix they have nothing to avoid: plain Cholesky; Cholesky taking the largest remaining diagonal entry as pivot, the pivoted form that reveals rank on a semidefinite matrix; and the two symmetric indefinite rules of the earlier essay, Bunch–Kaufman and the bounded rule. On a definite matrix the last two take no interchanges and no two-by-two blocks — their pivots are plain Cholesky’s to the last bit on all ninety matrices — so they drop out of the comparison, which is itself the first answer: neither rule has anything to offer a definite matrix that Cholesky does not.

No pivot of any factorisation is below the smallest eigenvalue, on any matrix. And none is a usable shift: unpivoted, the smallest pivot is never nearer the smallest eigenvalue than the second, on any of the ninety. Taking the largest diagonal first brings it nearer — on two matrices of the thirty in the most ill-conditioned family it crosses the midpoint — but on the other eighty-eight it does not.

Six to ten pushes, not four #

The search is the earlier essay’s, mirrored, and it rests on the fact that an eigenvalue count that cannot be slightly wrong measured: the signs of an elimination’s pivots count eigenvalues below a shift exactly, so whether a shift is below the spectrum is a question with a yes-or-no answer and no rounding in it. Start at the smallest pivot, and shrink the shift by a factor of 1.25 until the factorisation of $A - \sigma I$ succeeds — all pivots positive, which for a definite matrix’s shift means $\sigma$ is below $\lambda_{\min}$. The prediction said a median of four.

The medians are nine on the shifted family, six on the Gram matrices with forty-eight rows and ten on those with thirty. The fewest on any matrix is four and the most fourteen. Every step is one more factorisation of a 24 × 24 matrix, so the search that the earlier essay found cheap on indefinite matrices — two or three factorisations — costs three times as many on definite ones.

The number of steps is just the logarithm of the pivot’s distance above the spectrum, in units of $\log 1.25$, and that distance is large. The smallest pivot is 4.1 to 11.3 times the smallest eigenvalue on the shifted family, with a median of 6.4; 2.1 to 5.5 times on the forty-eight-row Gram matrices; and 3.9 to 21.6 times on the thirty-row ones. The prediction’s four steps corresponds to a factor of 2.4, which only the best-conditioned matrices approach.

The pivot follows the harmonic mean #

The prediction’s reason was that “a definite matrix’s pivots sit near its diagonal while its smallest eigenvalue sits far below it”. The first part is not quite what happens.

The pivots start near the diagonal — the first pivot is the first diagonal entry — and then fall, by a factor of 2.7 on the matrix the dial opens on. On the Gram matrix the dial opens on, the diagonal entries are all near one; the unpivoted pivots drift down from about one to 0.37 at the twenty-third position, and stop there. They do not fall to the smallest eigenvalue, 0.089, or even to the second, 0.148. They stop at the harmonic mean of the spectrum, 0.45.

That is an identity once it is looked at the right way. The last pivot of a Cholesky factorisation is one over the last diagonal entry of $A^{-1}$, and writing $A^{-1}$ in the eigenvector basis gives

$$d_n = \Big(\sum_i \frac{v_{ni}^2}{\lambda_i}\Big)^{-1},$$

where $v_{ni}$ is the last component of the $i$-th unit eigenvector. That is a weighted harmonic mean of the eigenvalues, with weights $v_{ni}^2$ that sum to one. The same holds for every pivot, with the eigenvalues and eigenvectors of the trailing block it is taken from. If the weights were all equal it would be the spectrum’s harmonic mean exactly; for a matrix whose eigenvectors are not aligned with the coordinates they are roughly equal, and the pivot is roughly the harmonic mean. It reaches $\lambda_{\min}$ only if the last coordinate vector is the smallest eigenvector, which for a matrix with random eigenvectors is as unlikely as anything gets.

Across all ninety matrices the smallest pivot is between 0.45 and 1.44 times the harmonic mean, while it is between 2.1 and 21.6 times the smallest eigenvalue. The harmonic mean of a spectrum is dominated by its small eigenvalues, so it sits much nearer the bottom than the arithmetic mean does — but not at the bottom. On the shifted family it is a median of 8.9 times $\lambda_{\min}$, because the shifted semicircle has only a thin tail of eigenvalues near its left edge; on the forty-eight-row Gram matrices 4.4 times; on the thirty-row ones 9.4 times. How far the smallest pivot sits above the spectrum is a statement about how thin the spectrum’s bottom tail is, and the harmonic mean is the number that measures it.

Taking the largest diagonal first pushes the last pivots lower, below the harmonic mean — 0.2 to 0.8 times it — because choosing the largest remaining diagonal each time leaves for the end the coordinates the matrix weights least, which tend to be the ones whose eigenvector weights lean towards the small eigenvalues. It is the mechanism that makes pivoted Cholesky reveal rank on a semidefinite matrix, and on a definite matrix with a well-separated bottom it helps by a factor of about two: medians of 3.9, 2.7 and 3.7 times $\lambda_{\min}$ instead of 6.4, 3.7 and 8.1. That is still not a usable shift.

Zero was always below the spectrum #

The pushing search was the right device for the indefinite case for a reason that does not apply here. An indefinite matrix’s most negative eigenvalue can be anywhere below zero, and the search has no floor to aim at: it can only walk outwards from the pivot until it passes the spectrum’s end. A definite matrix comes with a floor. Zero is below every eigenvalue, so the interval from zero to the smallest pivot is a bracket for $\lambda_{\min}$ before anything is computed.

Halving towards zero — $\sigma = d/2, d/4, \ldots$ until the factorisation succeeds — uses it. The medians fall to three, two and three and a half, never more than five; a little more than a third of the pushes. With the interval from the last failed shift to the first successful one as a bracket, a bisection step or two then places the shift as close to $\lambda_{\min}$ as the earlier essay’s bisection did.

The halving search is not free, and the figure says what it costs. The rate at which inverse iteration converges depends on where the shift sits relative to the two smallest eigenvalues, and the gap decides the eigenvector is the reminder that the gap between those two decides more than the rate: it decides how well the vector the iteration converges to is determined at all. A shift that lands well below the spectrum converges more slowly than one just below it, because inverse iteration’s rate is the ratio of the shift’s distances to the smallest and second eigenvalues. Halving lands further below: inverse iteration then needs a median of 5, 10 and 7 steps on the three families against 4, 8 and 5 after the finer shrinking. Each of those steps is a pair of triangular solves with the factor already in hand, about $2n^2$ operations against a factorisation’s $n^3/3$ — a quarter of the cost at this order, and a smaller share the larger the matrix — so trading four to six factorisations for one to three solves is a clear gain. On the forty-eight-row Gram family, whose gap between the two smallest eigenvalues is the narrowest — a median of half the smallest eigenvalue — the slowest matrix needs 65 steps after halving, which is the reason to spend one bisection step before starting.

A factorisation with nothing to pivot for is the account of why Cholesky needs no pivoting for stability on a definite matrix. This is the other side of the same structure: with no pivoting decisions to make, the factorisation also has nothing to say about the bottom of the spectrum beyond the harmonic mean, and the search for the bottom has to be done with inertia. The shift that stops at the first right count met the limit of that from the other side, on a saddle-point solver whose shift stops wherever the inertia first reads right: a count answers whether a shift is past the spectrum exactly, and says nothing about how far past.

What the shift is for #

None of this makes the smallest pivot useless; it makes it a different number from the one wanted. As an upper bound on the bottom of the spectrum it is free and guaranteed, which is more than can be said of anything a Krylov method returns in its first few steps, and the earlier essay on Cholesky, a curvature direction the factors cannot refine, found how little an indefinite factorisation’s own directions are worth once they are refined. What the definite pivot gives is the top of a bracket whose bottom is zero, and that is the right object for a search whose currency is factorisations: three halvings of a bracket are three factorisations, and the bracket is there before the first.

The practical rule that falls out is short. On a definite matrix, take the smallest pivot of the factorisation already computed, halve it until the shifted factorisation succeeds — each attempt a Cholesky factorisation that either completes or meets a nonpositive pivot — and spend one bisection step between the last failure and the first success before starting inverse iteration. On these ninety matrices the halving alone is a median of two to three and a half factorisations against six to ten by shrinking a quarter at a time, and inverse iteration from the halved shift a median of five to ten solves; the bisection step’s saving in solves is the one number here not measured.

What the prediction got right #

The prediction’s mechanism — pivots near the diagonal, the smallest eigenvalue far below — gets the direction right and the size wrong. The pivots do start near the diagonal, but they do not stay there; they fall by a factor of two to three over the factorisation, to the harmonic mean, and the factor that remains is the spectrum’s own: harmonic mean over smallest eigenvalue. A prediction built on the diagonal would have overshot as far as the four pushes undershot, since the diagonal entries average to the arithmetic mean of the spectrum, which on these families is well above the harmonic one; a prediction built on the harmonic mean says six to ten, and that is what was measured.

The earlier essay’s indefinite matrices are the instructive contrast. There the most negative pivot was nearest the most negative eigenvalue on twelve matrices of thirty under one rule and five under the other, and sometimes below it — more negative than any eigenvalue — because an indefinite pivot is not a weighted mean of anything: the weights in the identity above are the same, but with eigenvalues of both signs the reciprocals can cancel, and a weighted harmonic mean of mixed-sign numbers can be anywhere. A definite factorisation’s pivots cannot do that, which is the reason they can never be below the spectrum, and the same reason they can never be far from its harmonic mean. Where the multipliers go found the indefinite rules’ size problems in L; on a definite matrix L is tame and D is predictable, and predictably wrong as a shift.

What ninety matrices do not show #

The families are all dense with eigenvectors spread over every coordinate, which is what makes the weights nearly equal. A matrix whose smallest eigenvector is concentrated on a few coordinates — a discretised operator with a localised mode, a graph Laplacian with a weakly connected vertex — would put a large weight on $\lambda_{\min}$ for a pivot at one of those coordinates, and if that coordinate is eliminated last the smallest pivot would be close to $\lambda_{\min}$. Pivoting by the largest diagonal tends to do exactly that, and on such matrices it may be a usable shift. The order here is 24; the identity holds at any order, but the harmonic mean’s distance above the bottom grows with the order for most families, and so would the pushes.

Still open: a localised bottom, and a shift from the last pivots together #

A localised smallest eigenvector. The prediction with a sign is that on a path-graph Laplacian with one weak edge, shifted to be definite, the smallest pivot under largest-diagonal-first pivoting lands within a factor of 1.2 of $\lambda_{\min}$ and is a usable shift unpushed, while the unpivoted pivot stays near the harmonic mean — the first family on which the pivoting rule decides whether the factorisation knows the bottom of the spectrum.

The last few pivots as an estimate. Each of the last $k$ pivots is a harmonic mean with its own weights. The inverses of the last $k$ pivots are diagonal entries of inverses of nested trailing blocks, and their differences are what a Lanczos run on $A^{-1}$ would read off as moments. The prediction is that the smallest root of the $2 \times 2$ moment matrix built from the last three pivots is within a factor of 1.5 of $\lambda_{\min}$ on four matrices in five of the thirty-row Gram family — which would make the factorisation’s last pivots a two-step shift for no extra factorisation.