Neither sparse nor dense

The accuracy worth paying for

Used as a preconditioner, a hierarchical representation gets better at every accuracy — the iteration count falls monotonically all the way to the tightest tolerance. The total work does not. Its minimum sits at a rank-one preconditioner on an easy problem and six decades further along on a hard one.

Worth reading first: The exact answer to a nearby problem · Changing the condition number on purpose · Which pairs are allowed to be small.

The previous essay established that compressing a matrix is perturbing the problem. A perturbed problem is what a preconditioner is, so the obvious thing to do with a deliberately bad representation is to put it inside an iteration and let the iteration repair it.

That works. What it costs is the interesting part.

What a hierarchical preconditioner costs, against how accurate it is, at κ = 2.09e+1The iteration is conjugate gradients on the true matrix and the preconditioner is a solve with the representation, so what is measured is how much of the compression the outer loop repairs and what the repair costs. The count falls monotonically — 13, 8, 6, 3, 2, 2 steps against 41 with no preconditioner — and it keeps falling all the way to the tightest accuracy, which is the number a solver prints and the one everybody watches. The total work, one dense product with A plus one solve with the compressed matrix per step, is the other curve, and its minimum is at ε = 0.5 where the preconditioner has rank 1. Past that point every extra digit in the preconditioner is bought and not needed: at the tightest setting the run costs 2.42 times the cheapest, having taken 11 fewer steps to get there.-9-7-5-3-1110⁵10⁶10⁷log₁₀ ε of the preconditionermultiplications for the whole solveno preconditioner: 41 stepsiteration count, rescaledmultiplications spenttwo curves, opposite waysκ21steps, no preconditioner41cheapest ε0.5its rank1tightest ⁄ cheapest2.4the count is what is printedand the cost is what is spent
Fig. 1 Two curves going opposite ways. The dashed one is the iteration count and the solid one is the work, against how accurate the preconditioner is.

The set-up

Conjugate gradients on the true matrix A, preconditioned by a solve with AH. The iteration never sees the representation as the operator — every step forms one product with the real A — so the answer it converges to is the answer to the real problem, at whatever tolerance the iteration is asked for, and the compression cannot affect it — the residual it stops on is the real one.

What the compression affects is the rate. AH is a nearby matrix, so M⁻¹A is near the identity, and how near is exactly the backward error the previous essay measures.

At 128 unknowns and κ = 20.9, with 41 steps unpreconditioned:

ε of the preconditioner rank steps multiplications
5·10⁻¹ 1 13 4.81·10⁵
10⁻¹ 2 8 5.23·10⁵
10⁻² 3 6 6.26·10⁵
10⁻⁴ 6 3 7.13·10⁵
10⁻⁶ 8 2 7.43·10⁵
10⁻⁸ 10 2 1.16·10⁶
none — 41 6.72·10⁵

The step column falls monotonically. The work column rises monotonically. And the best entry in the work column is the rank-one preconditioner, which is a representation nobody would accept as an answer to anything.

Why the two columns disagree

Each step of the preconditioned iteration costs one product with A, which is n² and fixed, plus one solve with AH, which is not fixed. The solve’s cost grows with the ranks, and the ranks grow with the accuracy at a rate that field measures.

So the total is

steps × (n² + cost of one solve with AH) ,

and the accuracy pushes the first factor down and the second up. There is no reason for the product to be monotone in either direction, and it is not.

At κ = 20.9 the first factor runs out of room almost immediately. Thirteen steps is already cheap; getting to two saves eleven products with A, and buying them costs a preconditioner of rank ten where rank one would have done, whose solve is several times more expensive. The trade stops paying after the first step of it.

The iteration count is the number a solver prints and the number every convergence plot in this collection is drawn against, and here it is the wrong number. That is not a criticism of iteration counts — they are the right currency when the preconditioner is fixed, which is the situation every other preconditioning essay on this site is in. It is a warning about what happens when the preconditioner acquires a knob — the same warning the sequence field’s rebuild period carries about a rule written in iterations when the decision is half about the code.

The optimum moves, and with exactly one thing

If the trade is between iteration count and cost per iteration, the balance ought to shift when the iterations get more expensive to avoid — which is to say when the problem gets harder.

The matrix has a shift in it that was believed to do nothing to the geometry. Every distance between every pair of clusters is unchanged, so every rank in the partition is unchanged, and the condition number moves over three decades:

shift κ steps, unpreconditioned cheapest ε its rank its cost
1/n 2.09·10¹ 41 5·10⁻¹ 1 4.81·10⁵
4/n 1.75·10² 72 5·10⁻¹ 1 5.92·10⁵
16/n 1.53·10³ 105 10⁻⁶ 6 4.27·10⁵
64/n 1.07·10⁴ 119 10⁻⁸ 5 3.45·10⁵

Six decades of tolerance between the easiest problem and the hardest. On the easy one the best preconditioner in the comparison is rank one; on the hard one it is the tightest accuracy on the sweep.

What a hierarchical preconditioner costs, against how accurate it is, at κ = 1.53e+3The iteration is conjugate gradients on the true matrix and the preconditioner is a solve with the representation, so what is measured is how much of the compression the outer loop repairs and what the repair costs. The count falls monotonically — 17, 17, 8, 4, 2, 2 steps against 105 with no preconditioner — and it keeps falling all the way to the tightest accuracy, which is the number a solver prints and the one everybody watches. The total work, one dense product with A plus one solve with the compressed matrix per step, is the other curve, and its minimum is at ε = 10⁻⁶ where the preconditioner has rank 6. Past that point every extra digit in the preconditioner is bought and not needed: at the tightest setting the run costs 1.34 times the cheapest, having taken 0 fewer steps to get there.-9-7-5-3-1110⁵10⁶10⁷log₁₀ ε of the preconditionermultiplications for the whole solveno preconditioner: 105 stepsiteration count, rescaledmultiplications spenttwo curves, opposite waysκ1527steps, no preconditioner105cheapest ε10⁻⁶its rank6tightest ⁄ cheapest1.3the count is what is printedand the cost is what is spent
Fig. 2 The third row drawn. At κ = 1.53·10³ the step counts across the accuracy sweep are 17, 17, 8, 4, 2 and 2 against 105 unpreconditioned, and the totals are 6.29·10⁵, 6.29·10⁵, 5.23·10⁵, 4.77·10⁵, 4.27·10⁵ and 5.72·10⁵ — cheapest at ε = 10⁻⁶ with rank 6.

The total curve has a genuine interior minimum, and it is shallow. The six totals span 6.29·10⁵ down to 4.27·10⁵ and back up to 5.72·10⁵ — the worst is 1.47 times the best across six decades of tolerance, with the optimum 32% below the two loosest settings and 25% below the tightest. So the optimum is real, it is interior rather than at an end, and choosing it badly costs less than a factor of one and a half.

That is worth saying beside the table’s headline, because the table invites the opposite reading. Six decades of movement in the argmin sounds like a parameter that must be tuned per problem; a factor of 1.47 in the value says it is a parameter that mostly does not matter, provided it is somewhere in the middle. The expensive mistake is not picking the wrong ε — it is picking no preconditioner at all, where the step count is 105 rather than 2.

That was written here as a statement about conditioning alone, on the reasoning that the shift is not a geometric quantity and this field’s ranks are functions of geometry. The next section is why that is only half true, and what the other half does to the table.

The shift moves the ranks as well

The shift is added to the kernel’s own argument — the entry is f(|x − y| + δ) — so a larger δ moves every evaluation away from the singularity and makes the kernel smoother. Whether that matters depends on how δ compares to the distances the partition cares about, and at the far end of this sweep it does not compare favourably: δ = 64/n is 0.5 on a unit interval, which is half the domain rather than a shift on it.

Measured, the maximum rank the partition asks for at each tolerance:

at δ = 1/n the ranks are 1, 2, 3, 6, 8, 10 across ε from 5·10⁻¹ to 10⁻⁸. At 4/n: 1, 2, 3, 5, 7, 9. At 16/n: 1, 1, 2, 4, 6, 7. At 64/n: 1, 1, 2, 3, 4, 5.

At the tightest tolerance the rank halves across the sweep, and at 10⁻⁴ it runs 6, 5, 4, 3. So the ranks do move, and by a factor rather than by a rounding.

Which means the table varies two things

That changes what the previous section’s table is evidence for, and it is worth being exact rather than retracting it.

The optimum does move from ε = 5·10⁻¹ to ε = 10⁻⁸ across the sweep, and it moves for two reasons rather than one. The problem gets harder, so there are more iterations to be bought and a tighter preconditioner has more to save — which is the mechanism the section describes and it is real. And a tight preconditioner gets cheaper: at δ = 64/n a preconditioner at 10⁻⁸ costs a rank of 5 where at δ = 1/n it would have cost 10, so the second factor in the product is smaller at the far end of the table than at the near end.

Both push the optimum the same way, which is why the trend is clean and why the confounding is invisible from the table. Separating them would need a sweep that moves κ without moving the kernel’s smoothness — a different shift, applied to the spectrum rather than to the argument, in the way the barrier field’s diagonal does — and this representation does not have one.

What survives without qualification is the finding the essay opens with: the step count and the work disagree, and the work has an interior optimum. That is measured at a single shift, in the first table, where nothing is varying but ε. The direction the optimum moves with conditioning survives too, because both mechanisms agree on it. What does not survive is the attribution — the size of the move is not a measurement of how much harder the problem got.

That is the ordinary hazard of a knob that was believed to be orthogonal, and the site has a name for the check that catches it: vary the thing that was supposed not to move. Here it is one call to the same sweep asking for ranks instead of costs, and it would have been asked for if the claim had been stated as a measurement rather than as a reason.

The accuracy worth paying for in a preconditioner, against the condition number of the problemFour complete cost sweeps, each reduced to the ε that minimised the total work, against the κ of the problem it was run on. The shift that moves κ over 2.7 decades leaves every distance between every pair of clusters exactly where it was — and the ranks the partition asks for at a fixed tolerance do not stay put: they fall 10, 9, 7, 5 across the sweep, which is about half. So the shift is not the purely spectral knob it was chosen to be; adding it to the diagonal makes the kernel decay faster as well, and the second effect is doing part of the work the first is credited with. The optimum walks 0.5 → 0.5 → 10⁻⁶ → 10⁻⁸ — six decades of accuracy, bought because the problem got harder and for no other reason. On the easiest problem the best hierarchical preconditioner in the comparison has rank one; on the hardest it is the tightest one on the sweep. How accurate an approximate inverse should be is a question with an answer, and the answer is not in the matrix's structure.10¹10²10³10⁴10⁵-9-7-5-3-11condition number of the problemlog₁₀ ε at the cheapest totalrank 1rank 1rank 6rank 5the knob answers to κκ = 210.5κ = 1.1·10⁴10⁻⁸rank at the easy end1rank at the hard end5ranks at 10⁻⁸, spread5the geometry did not moveand the right accuracy did
Fig. 3 Four complete sweeps reduced to four points: the accuracy that minimised the work, against the condition number of the problem it was run on.

What the rank-one preconditioner is

A representation at ε = ½ deserves a paragraph, because rank one sounds like a joke and the object is not one.

Every admissible block is replaced by its single dominant term: the outer product of the block’s leading left and right singular vectors, scaled by the leading singular value. Every dense block is kept exactly. So the preconditioner is the exact near field plus a rank-one summary of every far-field interaction, and it reproduces the matrix to 7.8 per cent.

Physically that is the most natural approximation there is. The interaction between two well-separated clusters, to leading order, is the total of one cluster acting on the total of the other — a monopole against a monopole — and rank one is exactly that. The second term is the dipole, the third the quadrupole, and the columns-per-decade measurement in this field’s second essay is the statement that each multipole order buys a couple of digits.

So the first row of the table is not a degenerate corner of the sweep. It is the classical approximation the whole subject started from, and on an easy problem it is the best preconditioner in the comparison. What the later rows buy is multipole orders nobody needed.

Two consequences worth carrying. The first is that a hierarchical preconditioner at low accuracy is cheap to build as well as to apply: rank one per block means one dominant singular pair per block, which a few steps of a power iteration find without a decomposition. The second is that the whole first half of this field’s storage story — the partition, the tree, the admissibility test — is doing its work at ε = ½ exactly as much as at 10⁻⁸. The structure is what is expensive to decide and it is decided once.

Which is the same shape as four other knobs here

This collection has a field about what happens when a method acquires a parameter and it found that four of the site’s knobs are one knob: a regularisation parameter, a truncation rank, a stopping index and a filter width all control how much of the data to believe, and how much is worth believing belongs to the problem rather than to the method.

This is a fifth of that kind and it belongs in the same sentence with one word changed. How accurate an approximate inverse should be is a question whose answer is in the problem — specifically in κ — and not in the approximate inverse.

What is different, and what makes it worth a separate essay, is that the four knobs in that field all control an accuracy of the answer and this one does not. The answer’s accuracy is set by the iteration’s stopping tolerance and is the same at every row of both tables above. What the knob controls is only how much work the iteration does to get there, so the trade is entirely inside the cost and none of it is inside the result.

That is a rare and comfortable position. Every other parameter on this site has to be chosen against a quantity nobody can measure; this one has to be chosen against a quantity — total work — that the code can measure exactly, on the run it is doing, and adjust.

Against no preconditioner at all, which is the other comparison

Every row of both tables beats the unpreconditioned run, and by how much is worth reading because it is the part that makes the whole exercise worth doing.

At κ = 20.9 the unpreconditioned iteration takes 41 steps and 6.72·10⁵ multiplications; the best preconditioned run takes 13 and 4.81·10⁵. That is a saving of 28 per cent, which is real and unexciting. At κ = 1.07·10⁴ the unpreconditioned run takes 119 steps and 1.95·10⁶; the best preconditioned run takes 2 and 3.45·10⁵. That is a factor of 5.7.

So the preconditioner earns its keep in proportion to how hard the problem is, which is what a preconditioner is supposed to do and is worth confirming rather than assuming. It also means the two tables are answering two different questions and both are worth asking: is a hierarchical preconditioner worth having has the answer yes and increasingly so, and how accurate should it be has the answer that moves.

There is a third question the tables do not answer and it is the one the next essay takes. At every row above, the representation is being used as a preconditioner — applied inside an iteration that converges to the true answer. It could instead be used as a solver, applied once, returning an answer at the backward error the previous essay measured. On this matrix at 128 unknowns that costs 5.66·10⁵ multiplications, which is worse than four of the six rows in the first table, and the comparison flips completely by 512 unknowns.

What a code should actually do

The measurement suggests an obvious policy and the obvious policy is nearly right.

Since the total is measurable and the curve has one minimum, a code can start with a cheap preconditioner, solve, and refine the representation only if the iteration count came out high. The refinement is cheap in a way worth noticing: a hierarchical representation compressed at 10⁻² can be re-compressed at 10⁻⁴ from the blocks it already has, if the original block decompositions were kept, so a code that is prepared to spend the memory can walk down the curve without recomputing.

Two things spoil this and both are worth stating.

The first is that the setup cost is not on the table. Every row above prices only the iterations; the representation had to be built, and at a tight tolerance building it costs more. This collection has an essay in the sequence field about exactly that failure — a rule that reads a free measurement and beats every constant, except that the free measurement cannot see what a rebuild cost, and at the ends of the range the setup cost is the whole decision. The same hole is here.

The second is that the right preconditioner accuracy depends on how many right-hand sides there are. One solve amortises the setup over one iteration count; a hundred solves with the same matrix amortise it over a hundred, and then a tighter representation is worth building. The sequence field is the field about that, and this essay’s table is the one-right-hand-side column of a larger one.

The half of the sweep that is monotone, and the half that is not

There is a detail in the four condition-number sweeps worth pulling out, because it is the thing that makes the policy in the last section workable at all.

The step column is monotone at every condition number. At κ = 20.9 it goes 13, 8, 6, 3, 2, 2; at κ = 1.07·10⁴ it goes 15, 15, 14, 5, 3, 2. It never rises. That is not obvious in advance — a preconditioner is not required to improve when it gets more accurate, and this collection has a whole essay about a stale preconditioner that improves the problem while getting worse at it — but it is what happens, and it means the code never has to worry about overshooting the accuracy in the direction of the iteration count.

The cost column is not monotone at any condition number and its minimum is interior at two of the four. That is what makes the search a search rather than a ratchet.

And there is a shape in the step columns worth naming. At the two harder problems the count is flat for the first two or three tolerances — 15, 15, 14 — and then falls off a cliff. A preconditioner too crude for a hard problem does nothing at all, several decades of refinement do nothing at all, and then one further decade does everything. A code that walked the curve from the cheap end and stopped when refining stopped helping would stop at the first row and never find the minimum.

So the honest policy is not refine while it helps. It is refine past the point where it stops helping, at least once, which is an uncomfortable rule and is what the measurement says.

What a hierarchical preconditioner costs, against how accurate it is, at κ = 1.75e+2The iteration is conjugate gradients on the true matrix and the preconditioner is a solve with the representation, so what is measured is how much of the compression the outer loop repairs and what the repair costs. The count falls monotonically — 16, 14, 7, 4, 3, 2 steps against 72 with no preconditioner — and it keeps falling all the way to the tightest accuracy, which is the number a solver prints and the one everybody watches. The total work, one dense product with A plus one solve with the compressed matrix per step, is the other curve, and its minimum is at ε = 0.5 where the preconditioner has rank 1. Past that point every extra digit in the preconditioner is bought and not needed: at the tightest setting the run costs 1.59 times the cheapest, having taken 14 fewer steps to get there.-9-7-5-3-1110⁵10⁶10⁷log₁₀ ε of the preconditionermultiplications for the whole solveno preconditioner: 72 stepsiteration count, rescaledmultiplications spenttwo curves, opposite waysκ175steps, no preconditioner72cheapest ε0.5its rank1tightest ⁄ cheapest1.6the count is what is printedand the cost is what is spent
Fig. 4 The second sweep, where the count is 16, 14, 7, 4, 3, 2 and the plateau at the cheap end has almost disappeared.

The refusal

The claim under test is one that is stated in nearly every paper on preconditioning and is nearly always harmless: that a preconditioner which takes fewer iterations is a better preconditioner.

It is harmless when the preconditioners being compared cost the same to apply, which is the usual case. It is false as soon as one of them has an accuracy knob, and the refusal is fed the two ends of the easiest sweep: the tightest preconditioner takes 2 steps against the cheapest one’s 13, and costs 1.16 million multiplications against 0.48 million. The assertion that the tighter one costs less is handed those numbers and it fails.

The general form is one this collection keeps meeting. A currency is only a currency while the exchange rate is fixed, and an iteration count is a currency for work exactly as long as an iteration costs a constant. Two other essays here are about the same substitution going wrong — a flop count that stopped predicting the time once the memory hierarchy appeared, and a residual that stopped predicting the error once the condition number did.

What links here

Computed from the collection, not written here: the essays that point at this one.

Reads more easily once this is understood

Essays that name this one as worth reading first.

Shares its objects with

Essays that name at least two of the same things, and that neither author linked.

Named objects

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

Backward errorCondition numberConjugate gradientsFlop countHierarchical matrixParameter choicePreconditioningStopping criterionToleranceWoodbury identity