Numerical Linear Algebra
where the algebra stops being the arithmetic
Eigenvalues, singular values, rank

Restarting is a filter

A restart throws away the Ritz values it does not want and begins again from a new starting vector. Written in the eigenbasis, that vector's components have been multiplied by a polynomial with its roots at the discarded values — measured component by component, and agreeing with the polynomial to rounding.
17 min read 8 figures Measured, not assertedTwo routes to a numberIdentical algebra, different arithmetic

Worth reading first: An eigenvalue that arrives twice · The rate the condition number predicts.

The Lanczos essay ends with three costs, all of them growing with the number of steps. The basis is m vectors and they have to be kept. Orthogonality is lost geometrically over about ten consecutive steps and never recovered. And the method returns twenty-five extra copies of thirteen eigenvalues in eighty steps on a forty-by-forty matrix, each copy accurate to 1.9·10⁻⁸ and none of them a second eigenvalue.

Restarting is the standard answer to all three, and the standard description of it — run m steps, keep the k best pairs, start again — makes it sound like a housekeeping measure. It is not. The restart is where the method does its aiming, and the aiming is a filter in exactly the sense the regularisation field’s filter factors are one.

What 3 shifts do to each direction of the starting vectorAmplification of each eigen-component against its eigenvalue, on a logarithmic vertical axis. The points are measured — the component of the vector along each eigenvector, before and after — and the curve is |Π(λ − θ)| with the roots at the 3 discarded Ritz values. They agree to 3.9·10⁻¹¹. The 3 wanted directions are amplified by at least 1 and everything else by at most 0.802.1357910⁻⁵10⁻⁴10⁻³10⁻²10⁻¹1eigenvalue λamplificationkeptdiscardeda filter on the starting vectormeasured against the polynomial3.9·10⁻¹¹worst kept direction1best discarded direction0.8the roots are the discarded Ritz valuesand the vertical lines are where they sit
Fig. 1 What three shifts do to each direction of the starting vector. The points are measured — the component of the vector along each eigenvector, before and after — and the curve is |Π(λ − θ)| with its roots at the three discarded Ritz values, drawn as the vertical lines. They agree to 10⁻⁹ at every eigenvalue. Drag the number of shifts and watch the polynomial acquire roots.

The mechanism, written out

Run m = k + p steps of Lanczos from a starting vector v. That produces m Ritz values; keep the k wanted ones — the largest, here — and call the other p unwanted. Then form

v⁺  =  Π_j (A − θⱼ I) v  /  ‖·‖

with the θⱼ the unwanted values, and start again from v⁺.

Write v in A’s eigenbasis as Σ cᵢ xᵢ. Each factor (A − θI) multiplies the ith component by (λᵢ − θ), so the whole product multiplies it by Π_j(λᵢ − θⱼ) — a polynomial evaluated at λᵢ, and the same polynomial for every component. Directions whose eigenvalues sit near a discarded Ritz value are multiplied by something small; directions near the wanted values are multiplied by something large.

That is what the figure shows and it is checkable because the matrix was built from its own eigenvectors: withSpectrum takes a list of eigenvalues and a fixed random orthogonal V and returns VΛVᵀ. So the components can be measured directly rather than inferred, and the assertion is that every one of them is scaled by the polynomial times a single shared normalisation — checked at all twenty eigenvalues, not at the extremes.

What the implicit method does instead, and why the difference is not visible here

The version above forms the product explicitly: p products with A per cycle, and a starting vector that has been multiplied by a polynomial in A.

The implicitly restarted method reaches the same vector without forming any of it. It applies p steps of shifted QR to the m×m tridiagonal, which — by the implicit Q theorem — performs the same filtering on the Krylov basis at a cost that does not involve A at all. That is cheaper and it is the same vector.

Drawing the explicit version is a deliberate choice, and this site has made it before. The Francis double shift draws the implicit step and then checks it against the explicitly formed (A − μI)(A − μ̄I) = QR, entry by entry, agreeing in absolute value to 2.2·10⁻¹⁵. The explicit form is the argument; the implicit form is an implementation of it, and separating the two is what lets a reader see what the machinery is for before meeting the machinery.

What the filter buys

12 restarts keeping 4 of 8, on a 40×40 matrixThe residual bound of the worst wanted eigenvalue and its true error, against the number of products with A. The bound falls from 1.21 to 1.08·10⁻¹³ across 12 cycles and 140 products, and the true error reaches 7.11·10⁻¹⁵. The basis is 8 vectors at every cycle and never grows.0183654729010812610⁻¹⁶10⁻¹³10⁻¹⁰10⁻⁷10⁻⁴10⁻¹products with Asizeresidual boundtrue errorbounded memorybasis vectors kept8products with A140worst error in the k wanted7.1·10⁻¹⁵the bound is free and the error is notand the basis never grows
Fig. 2 Twelve restarts keeping four Ritz values out of eight, on a 40×40 matrix. The residual bound falls from 1.21 to 1.1·10⁻¹³ across 144 products with A, and the true error of the worst wanted eigenvalue reaches 7.1·10⁻¹⁵. The basis is eight vectors at every cycle and never grows.

Eight vectors. An unrestarted run reaching the same accuracy on this matrix keeps every vector it has generated, and the whole reason the method exists is that at n = 10⁶ the difference between eight vectors and eighty is the difference between a method that runs and one that does not.

The comparison is drawn in products with A rather than in cycles, because that is the currency both methods spend: a cycle of m = 8 steps costs eight products plus p = 4 for the filter, so twelve cycles is 144 products. An unrestarted run of 144 steps on a 40×40 matrix would exhaust the space five times over — which is the other half of why restarting is not merely a memory measure.

Two properties of the convergence are worth separating.

The residual bound is free. βₘ|s_{m,k}| needs the last entry of the tridiagonal’s own eigenvector and no product with A at all. It is what a real code stops on, and it is what is plotted; the true error is drawn beside it only because this matrix was constructed and the answer is known.

And the progress is not monotone at every cycle, which is why the assertion is stated as a ratio across the whole run rather than step by step. A restart is not required to improve the worst wanted Ritz pair at every cycle, and claiming it does would be claiming something false.

The ghosts, as a side effect

lanczos.js measures what the plain three-term recurrence does over eighty steps: orthogonality lost by twelve orders of magnitude, and twenty-five duplicate eigenvalues that are accurate to 1.9·10⁻⁸ and are not there.

A restarted method with full reorthogonalisation inside each cycle cannot do that, and the reason is arithmetic rather than clever: the basis is eight vectors, reorthogonalising against eight vectors costs eight projections a step, and orthogonality lost geometrically over ten steps has no ten steps to be lost over. The ghosts are a consequence of a long recurrence, and there is no long recurrence.

Copies of each eigenvalue after 80 steps on a 40×40 matrix with a simple spectrumA bar per eigenvalue that came back more than once, showing how many times. The matrix has 40 distinct eigenvalues by construction; the run returned 25 extra copies of 13 of them, the most-copied arriving 5 times. Every copy is accurate to 1.9·10⁻⁸ relative, which is why nothing but the true spectrum could detect them.λ = 105 timesλ = 9.55 timesλ = 95 timesλ = 8.55 timesλ = 2.952 timesλ = 2.92 timeseigenvalues that arrived more than once — the matrix has 40 distinct onesa spectrum with the wrong multiplicitiesextra copies, no reorthogonalisation25extra copies, full reorthogonalisation0worst relative error among the copies1.9·10⁻⁸steps taken of 80 asked for, full40no arithmetic error was madeevery one of these is right to eight digits
Fig. 3 What the long recurrence produces, from the field’s own essay: extra copies of eigenvalues it has already found, each accurate to 1.9·10⁻⁸, and more of them the longer the run. A restarted method bounds the run and the ghosts with the same decision.

So restarting solves three problems with one mechanism, and only one of the three is what it is usually introduced for. That is worth noticing because it changes what the parameters are chosen for: k and p are usually described as a memory budget, and they are also the degree of the filter and the length of the recurrence.

Choosing the shifts, which is choosing the polynomial

The shifts here are the unwanted Ritz values, which is the standard choice and is called exact shifts. Read as a filter, it says: put the polynomial’s roots where the method’s own current estimate says the uninteresting part of the spectrum is.

That is a good choice and it is not the only one, and the filter reading is what makes the alternatives legible rather than arbitrary. Shifts placed on a line segment covering the unwanted region — Leja points, Chebyshev points — are choosing a polynomial small on an interval rather than at m points, and they behave differently when the unwanted spectrum is a continuum rather than a few isolated values. None of that is measured here; what is measured is that the mechanism is a polynomial, at which point the design space is the space of polynomials and not a list of tricks.

The measurement in the figure says how well the exact-shift choice does its job on this matrix: the worst wanted direction is amplified by 1.005 and the best discarded one by 0.802, a ratio of 1.25 in one cycle. That is not a dramatic separation — and it does not need to be, because the cycle repeats and the ratio compounds. Twelve cycles of 1.25 is a factor of fourteen, and the bound in the figure above falls by thirteen orders of magnitude over those twelve cycles, because the Ritz values sharpen as the vector improves and each cycle’s polynomial is better aimed than the last one’s.

The connection this essay exists to make

Three filters have now been measured on this site and they act on three different objects.

the method the filter acts on the parameter
Tikhonov, truncation the solution’s spectral components λ, K
conjugate gradients the solution, through the residual polynomial the step count
a restart the starting vector the shifts

The first two are the previous field’s subject and they suppress components of an answer. This one suppresses components of an input, and what it is aiming at is not accuracy but attention — which directions of the space the method will spend its next m steps on.

They are the same mathematics in the sense that matters: a polynomial in A applied to a vector, with the design question being where to put its roots. The site’s identical-algebra-different-arithmetic thread is usually about two computations of one quantity; this is one computation appearing in two fields with different names, which is the shape a synthesis phase is for.

The filter 8 conjugate gradient steps apply, measured and predictedFilter factors against the singular-value index. The factors measured off the iterate and the polynomial predicted from the recurrence coefficients agree to 6.9·10⁻¹⁴ and are drawn as one curve. It rises above one — 1.070 at its largest — and changes direction several times. Tikhonov's filter at the matching cutoff is monotone and never exceeds one.081624324048566400.250.50.7511.25index kfilter factoronezeroCG, both routesTikhonovone of these is not a weightlargest CG factor1.1measured vs predicted6.9·10⁻¹⁴largest Tikhonov factor0.97two routes to the same curveand a curve that goes above one
Fig. 4 The other polynomial, from the regularisation field. The conjugate gradient filter after eight steps, with its roots at that method’s Ritz values — the same construction, applied to the solution rather than to the starting vector, and drawn against the smooth filter it competes with.

What a restart cannot do

One limitation belongs here rather than in the what is left, because it is the subject of the next essay and because it is a consequence of the filter reading rather than an unrelated fact.

Everything above multiplies the components of a single starting vector. Whatever the polynomial does, the space the next cycle explores is the Krylov space of one vector — and a repeated eigenvalue’s eigenspace contributes exactly one direction to that, no matter what v is. So a restart can aim the method at a part of the spectrum and cannot make it see a multiplicity.

The refusal this essay’s family runs makes the same point from the other side: fed a starting vector whose component along the largest eigenvalue has been removed entirely, the claim that a Krylov method finds that eigenvalue anyway is required to fail. It does. A filter can suppress a direction, and it cannot create one.

What each reorthogonalisation costs and what it leaves, n = 40Two quantities for three variants on one logarithmic axis: the number of projections performed, and the departure from orthogonality left behind. Doing nothing costs nothing and leaves a basis that is not one. Full reorthogonalisation costs 820 projections and leaves 1.9·10⁻¹⁵. Selective costs 137 — 6.0 times fewer — and leaves 10·10⁻⁸, which is √u and is enough.010⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²10¹count, and ‖QᵀQ − I‖√u — semi-orthogonalitynoneselectivefullprojections‖QᵀQ − I‖two counters, three methodsghosts, none4ghosts, selective0selective is cheaper by6the cheap repair leaves a basis that is not orthogonaland that is the claim, not a shortfall
Fig. 5 The alternative to bounding the recurrence: reorthogonalise it. Selective reorthogonalisation is six times cheaper than full and leaves a basis at √u; a restart bounds the basis instead, and the two are the field’s two answers to the same problem.
Driving the subdiagonal to zero, with λ₄/λ₃ = 0.90A semi-logarithmic plot of the magnitude of the subdiagonal entry against iteration count for three shift strategies. The unshifted curve is a straight line; the two shifted curves plunge to the bottom of the plot within a few steps.0357010514017521024528010⁻¹⁶10⁻¹³10⁻¹⁰10⁻⁷10⁻⁴10⁻¹iteration|subdiagonal entry|no shiftRayleighWilkinsontwo routes to one raterate, from the spectrum0.9rate, measured0.9iterations, none / Wilkinson45symmetric 4×4, spectrum 8, 4, 2, 1.8the dashed line is the prediction
Fig. 6 The shift in the method the libraries run, which is the same idea with one root: a polynomial in A applied to move the iteration’s attention, chosen from the method’s own current estimate of the spectrum.
Three filters at their own optima, on one problemFilter factors against the singular-value index. Truncation at K = 21 is a step from one to zero; Tikhonov at λ = 0.0215 is a smooth descent; the 20-step conjugate gradient filter rises above one and changes direction 5 times. The three reach errors within a few per cent of each other while differing by up to 0.99 in individual factors.081624324048566400.250.50.7511.25index kfilter factortruncation K = 21TikhonovCGLS, 20 stepsone floor, three weightingstruncation vs Tikhonov0.77truncation vs CGLS0.99reversals in the CGLS filter5the same error, to a few per centby three different weightings
Fig. 7 Three filters from the regularisation field, for scale. All of them weight the components of an answer; the restart above weights the components of an input, and the design question — where to put the polynomial’s roots — is the same one.
‖QᵀQ − I‖ at every step of a Lanczos run, n = 40The departure from orthogonality of the Lanczos basis, against the step, on a logarithmic vertical axis. It sits at the level of rounding for 13 steps and then climbs by a factor of about seventeen a step for ten steps running, twelve orders of magnitude, before saturating. The dashed line marks the step at which the first Ritz value converged, which is step 11.051015202530354010⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²10¹step‖QᵀQ − I‖first Ritz value converges√uno reorthogonalisationfullthe two events are one eventstep orthogonality crosses √u14step the first Ritz value converges11consecutive steps at ≥5× growth10a drift would grow like a square rootthis grows like a geometric series and then stops
Fig. 8 What a long recurrence does, and what a restart therefore prevents: orthogonality lost geometrically over about ten consecutive steps, arriving with the first converged Ritz value rather than at a fixed step count.

What the machinery had to be checked against

ritzPairs writes the Lanczos recurrence out a second time, because a restart has to start from a given vector and this site’s own lanczos starts from a seeded random one. A second body for a routine the site already has is exactly what kit_dup_check exists to be suspicious of, and the suspicion is right: a copy that drifts is how a measurement ends up describing a bug.

So it is held against lanczos.js on the one input where the two must agree — same matrix, same starting vector, same number of steps, full reorthogonalisation in both — and the Ritz values agree to 5.3·10⁻¹⁵. That check is cheap and it is load-bearing: without it the convergence in the figure above could be measuring the copy rather than the method.

The filter measurement needs a second piece of machinery for the same reason. The components of a vector along each eigenvector are only available because withSpectrum builds the matrix from a prescribed spectrum and a fixed random orthogonal V — the same move exact.js makes with the Hilbert inverse and fft.js with the Kac–Murdock–Szegő one. On a measured matrix the components would have to be estimated by an eigensolver, and comparing an eigensolver’s output against a polynomial would be comparing two approximations.

The residual bound, which is what a real code stops on

Every convergence claim in this essay is drawn against βₘ|s_{m,k}| rather than against the true error, and the reason is that the bound is free: it needs the last entry of the tridiagonal’s own eigenvector and no product with A at all. A code that has no exact answer to compare against — every code — stops on that number.

The true error is drawn beside it only because this matrix was constructed. What the pair shows is how much slack the bound has: at the last cycle the bound is 1.1·10⁻¹³ and the error is 7.1·10⁻¹⁵, so the quantity a real run reads is about fifteen times pessimistic there. That is the sort of factor a stopping tolerance has to be set with in mind, and it is the sort of number that is almost never reported beside a convergence curve.

Choosing k and p, which is choosing three things at once

The parameters are usually presented as a memory budget: keep k, work in k + p, and pick p as large as the machine allows. The filter reading says that the same two numbers are doing two more jobs.

p is the degree of the filter. Each shift is a root, so a larger p is a sharper separation between the wanted and unwanted parts of the spectrum per cycle — and costs p products with A per cycle to apply.

k + p is the length of the recurrence, which is what decides how much orthogonality is lost inside a cycle before it is thrown away. A restarted method with a very long cycle has the ghost problem back.

And k is what the answer is. Asking for four eigenvalues and keeping four is a different method from asking for four and keeping eight, because the four extra Ritz pairs carry information the next cycle’s filter is aimed with.

Nothing here measures the trade between them, and the natural experiment — the same total number of products with A, spent as many short cycles or few long ones — is a figure this field could have and does not.

What is left

Thick restarting, which keeps the wanted Ritz vectors explicitly rather than re-deriving them from a filtered starting vector, and is what most current codes do. It reaches the same subspace by a different arithmetic route, and the interesting comparison would be their behaviour under rounding rather than their exact-arithmetic equivalence.

Restarting for interior eigenvalues, where the exact-shift choice is much less obviously right — the unwanted spectrum surrounds the wanted values on both sides, and a polynomial small everywhere except a window in the middle is a harder object than one small at one end.

And the count of restarts as a parameter. Everything here fixes k, p and the number of cycles by hand. The step count became a regularisation parameter in another field this phase; whether the cycle count behaves the same way — with an interior optimum, and a rule for finding it — is a measurement nobody here has made.

Where this goes next

λ = 105 timesλ = 9.55 timesλ = 95 timesλ = 8.55 timesλ = 2.952 timesλ = 2.92 timeseigenvalues that arrived more than once — the matrix has 40 distinct onesa spectrum with the wrong multiplicitiesextra copies, no reorthogonalisation25extra copies, full reorthogonalisation0worst relative error among the copies1.9·10⁻⁸steps taken of 80 asked for, full40no arithmetic error was madeevery one of these is right to eight digits An eigenvalue that arrives twice A matrix with forty distinct eigenvalues, handed to Lanczos for eighty steps, returns twenty-five extra copies of thirteen of them — the largest arriving five times. Every copy is accurate to 1.9·10⁻⁸ relative. No arithmetic error was made, nothing overflowed, and a caller counting eigenvalues gets the wrong multiplicity from a computation in which no individual number is wrong. 015304560759010512010⁻²10⁻¹110¹steprelative sizeleast error: 20discrepancy stop: 7errorresidualthe knob is an integerleast error, at step20error there0.14error at step 1206the residual falls at every stepthe error turns and keeps rising A parameter that counts steps The regularisation field's knob is a positive real number chosen by one of three rules. The iterative field's is an integer nobody called a knob — where to stop. On the same problem the best step is 20 and the best λ is 0.025, and they reach 0.1426 and 0.1406. 12345610⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²10¹Ritz value, largest firstdistance from 10one vectora block of twoa space, not a ratecopies found, one vector1copies found, block of two2Krylov dimension, one vector16the second copy is not in the spaceat any number of steps An eigenvalue one vector cannot see A matrix with an exactly doubled eigenvalue at 10. Twelve Lanczos steps find it once; twenty-four find it once, on a Krylov space of dimension 23 in a 24-dimensional problem. A block of two vectors finds it twice. This is not slow convergence — the second copy is not in the space.

Named objects

A flat tag is an object no other essay names yet.

Filter factorsInvariant subspaceKrylov subspaceLanczos algorithmReorthogonalisationRestartingRitz valuesThree term recurrence