Neither sparse nor dense

A knob calibrated in residuals

A formatted Cholesky has two numbers in it and only one of them is an accuracy. Across twelve trees — three sizes by four leaf sizes — the leaf moves the truncation count from 0 to 258 and moves the ranks of the blocks not at all, while the residual follows the tolerance at slopes between 1.022 and 1.046 and sits at about a tenth of it throughout.

Worth reading first: The rounding that was not the problem · Where the format starts paying · A block nobody can call sparse.

A hierarchical factorisation is configured by two numbers, and a code that has to produce an answer to a stated accuracy has to know which of them to turn. One is a tolerance, applied to every admissible block when it is compressed and again after every truncation the arithmetic performs. The other is a leaf size: how few unknowns a block has to hold before it is stored densely instead of being split again.

Both of them change the answer. Neither of them is named in the units the caller cares about, which are the units of the residual. So the practical question is not what the two parameters mean but whether they separate — whether one of them can be set from the machine and the other from the accuracy required, with no term connecting them.

That question is answered here on a grid rather than on a line. Four leaf sizes, three matrix sizes, five tolerances each: sixty formatted Cholesky factorisations of the same kernel, every one of them checked against a dense reference. The answer is nearly the clean one, and the two places it is not clean are both visible only because the grid has a second axis.

The residual of a formatted Cholesky on 256 unknowns, against the accuracy its blocks were compressed atA leaf of 16 on 256 unknowns means 34 truncations inside the factorisation, and the ranks of its blocks run to 4, 7, 9, 12, 14 across the sweep. The residual tracks the tolerance at a slope of 1.027 over ten decades: the knob does what a knob should, and nothing in the approximate arithmetic bends it. The ratio to the representation's own error stays at 0.93, 1.83, 0.88, 0.91, 0.91 — one excursion above one, at 10⁻⁴, where the tolerance sits between two singular values of one block and the rank it takes is a rounding of a decision rather than a measurement.10⁻¹¹10⁻⁹10⁻⁷10⁻⁵10⁻³10⁻¹10⁻¹²10⁻¹⁰10⁻⁸10⁻⁶10⁻⁴10⁻²ε the blocks were compressed at‖A − LLᵀ‖ ⁄ ‖A‖the factorisationthe knob, on a factorisationtruncations34residual at 10⁻⁴2.2·10⁻⁵residual at 10⁻¹⁰9.9·10⁻¹²slope1worst ratio to the representation1.8ten decadesand a slope of one
Fig. 1 A formatted Cholesky on 256 unknowns with a leaf of 16. Thirty-four truncations happen inside it, the ranks its blocks take run 4, 7, 9, 12, 14 across the sweep, and the residual follows the tolerance at a slope of 1.027 over ten decades. The slider changes the leaf size.

The leaf moves the truncation count by a factor of forty-nine

The number of approximate operations inside the recursion is not something anybody sets. Every Schur update is a low-rank addition pushed down the tree, it reaches every admissible block below it, and each of those is truncated back. So the count is a consequence of the shape of the tree, and the leaf size is what fixes that shape.

At 256 unknowns with a leaf of 64 the tree is two levels deep and the factorisation performs two truncations. At a leaf of 8 it is five levels deep and performs ninety-eight. Same matrix, same tolerance, forty-nine times as many approximate operations — which is the reason the essay that first measured this cost was worth writing, because the received warning about formatted arithmetic is that a count like that has to show up somewhere.

The residual of a formatted Cholesky on 256 unknowns, against the accuracy its blocks were compressed atA leaf of 64 on 256 unknowns means 2 truncations inside the factorisation, and the ranks of its blocks run to 4, 7, 9, 12, 14 across the sweep. The residual tracks the tolerance at a slope of 1.022 over ten decades: the knob does what a knob should, and nothing in the approximate arithmetic bends it. The ratio to the representation's own error stays at 0.95, 0.95, 0.96, 0.93, 0.93 — one excursion above one, at 10⁻⁶, where the tolerance sits between two singular values of one block and the rank it takes is a rounding of a decision rather than a measurement.10⁻¹¹10⁻⁹10⁻⁷10⁻⁵10⁻³10⁻¹10⁻¹²10⁻¹⁰10⁻⁸10⁻⁶10⁻⁴10⁻²ε the blocks were compressed at‖A − LLᵀ‖ ⁄ ‖A‖the factorisationthe knob, on a factorisationtruncations2residual at 10⁻⁴5.2·10⁻⁶residual at 10⁻¹⁰9.5·10⁻¹²slope1worst ratio to the representation0.96ten decadesand a slope of one
Fig. 2 The shallowest tree at this size: a leaf of 64, two levels, two truncations. The residual runs 0.0014, 5.2·10⁻⁶, 9.7·10⁻⁸, 9.1·10⁻¹⁰ and 9.5·10⁻¹² across the five tolerances, at a slope of 1.022.

The count is worth stating exactly, because it turns out to be a function of one number rather than of two. Writing the depth of the tree as the number of levels between the root and a leaf, the twelve configurations measured here produce six distinct counts — 0, 2, 10, 34, 98 and 258 at depths of one through six — and each of those counts appears at every size that can reach that depth. Two hundred and fifty-six unknowns with a leaf of 16 and five hundred and twelve with a leaf of 32 are both four levels deep, and both perform exactly thirty-four truncations. One hundred and twenty-eight with a leaf of 8 makes it three configurations at the same count, across a factor of four in the size of the matrix.

The closed form the six counts satisfy is (L − 2)·2ᴸ + 2 at depth L, exactly, at every one of the six depths. It is a statement about the tree and not about the arithmetic: the count of admissible blocks touched by the Schur updates of a two-way recursion, and nothing in it knows what the entries are. That is worth having because it means the count can be predicted before anything is factorised, which is the opposite situation from an accumulation whose length is a parameter of the experiment.

The residual of a formatted Cholesky on 256 unknowns, against the accuracy its blocks were compressed atA leaf of 8 on 256 unknowns means 98 truncations inside the factorisation, and the ranks of its blocks run to 4, 7, 9, 12, 14 across the sweep. The residual tracks the tolerance at a slope of 1.043 over ten decades: the knob does what a knob should, and nothing in the approximate arithmetic bends it. The ratio to the representation's own error stays at 0.83, 1.78, 0.85, 0.81, 0.86 — one excursion above one, at 10⁻⁴, where the tolerance sits between two singular values of one block and the rank it takes is a rounding of a decision rather than a measurement.10⁻¹¹10⁻⁹10⁻⁷10⁻⁵10⁻³10⁻¹10⁻¹²10⁻¹⁰10⁻⁸10⁻⁶10⁻⁴10⁻²ε the blocks were compressed at‖A − LLᵀ‖ ⁄ ‖A‖the factorisationthe knob, on a factorisationtruncations98residual at 10⁻⁴2.2·10⁻⁵residual at 10⁻¹⁰10⁻¹¹slope1worst ratio to the representation1.8ten decadesand a slope of one
Fig. 3 The deepest tree at this size: a leaf of 8, five levels, ninety-eight truncations. The residual runs 0.0023, 2.2·10⁻⁵, 1.5·10⁻⁷, 1.1·10⁻⁹ and 1.0·10⁻¹¹, at a slope of 1.043 — the same line, lifted.

And it does not move the ranks at all

The second column of the figure is the one that decides the storage, and it is the surprise. Across all four leaf sizes at 256 unknowns the ranks the blocks take are 4, 7, 9, 12, 14 — the same five integers, one per tolerance, at a leaf of 8, of 16, of 32 and of 64. Forty-nine times as many truncations, and not one of them changes what any block keeps.

The mechanism is structural rather than numerical, and it is worth naming because it explains why the null result is exact rather than approximate. The largest rank anywhere in the tree is attained at the topmost admissible block, which is the coupling between the first half of the unknowns and the second half. That block exists at every leaf size, because the first split of the recursion happens before the leaf size is ever consulted; it is the same block, of the same dimensions, holding the same entries, compressed at the same tolerance. A smaller leaf adds blocks below it, and the blocks below it are smaller and have lower ranks. Nothing a deeper tree does can raise the maximum, and the measurement says nothing lowers it either.

So the leaf size is not merely accuracy-neutral in the residual, which was the earlier finding. It is neutral in the ranks, which is a stronger statement and a more useful one, because the ranks are what a storage estimate is built from. A code sizing its buffers from a sample factorisation at one leaf can carry that number to another leaf unchanged. What it cannot carry is the number of blocks, which is what the partition-against-storage trade is about, and which moves in the opposite direction from the dense fringe.

This is also the sense in which a rank here is not a knob at all. Rank is a decision when a threshold is laid across a decaying sequence, and the decision is made by the tolerance; the leaf size is not in that sentence anywhere.

The ranks belong to the problem, not to either knob

If the ranks were a property of the format they would move when the format moved, and they do not. The remaining possibility is that they are a property of the matrix, and the third axis of the grid is what tests it.

The residual of a formatted Cholesky on 128 unknowns, against the accuracy its blocks were compressed atA leaf of 16 on 128 unknowns means 10 truncations inside the factorisation, and the ranks of its blocks run to 3, 6, 8, 10, 12 across the sweep. The residual tracks the tolerance at a slope of 1.028 over ten decades: the knob does what a knob should, and nothing in the approximate arithmetic bends it. The ratio to the representation's own error stays at 0.97, 1.79, 0.90, 0.98, 0.97 — one excursion above one, at 10⁻⁴, where the tolerance sits between two singular values of one block and the rank it takes is a rounding of a decision rather than a measurement.10⁻¹¹10⁻⁹10⁻⁷10⁻⁵10⁻³10⁻¹10⁻¹²10⁻¹⁰10⁻⁸10⁻⁶10⁻⁴10⁻²ε the blocks were compressed at‖A − LLᵀ‖ ⁄ ‖A‖the factorisationthe knob, on a factorisationtruncations10residual at 10⁻⁴2.1·10⁻⁵residual at 10⁻¹⁰9.9·10⁻¹²slope1worst ratio to the representation1.8ten decadesand a slope of one
Fig. 4 The same leaf on half the unknowns: ten truncations rather than thirty-four, and the ranks fall to 3, 6, 8, 10, 12 across the sweep. The residual is unmoved — 0.0017, 2.1·10⁻⁵, 1.3·10⁻⁷, 1.0·10⁻⁹ and 9.9·10⁻¹² at a slope of 1.028.

Halving the matrix takes the ranks from 4, 7, 9, 12, 14 to 3, 6, 8, 10, 12. Doubling it takes them to 4, 7, 10, 13, 16. Those are three different columns of numbers for one parameter that neither knob controls, and they move in the direction the geometry says they should: the top-level block gets larger, its singular values decay at the same rate against a longer sequence, and the tolerance cuts it further along.

The rate is modest and it is not constant across the sweep. At the tightest tolerance measured the maximum rank runs 12, 14, 16 for a fourfold change in the number of unknowns — two columns per doubling. At the loosest it runs 3, 4, 4, which is one column across the whole range and then nothing. So the growth is logarithmic in the size and its constant depends on the tolerance, which is the same shape the exchange rate between digits and columns measures on a single block and at a lower rate, because the number here is a maximum over a tree of blocks rather than one block’s own count.

The residual of a formatted Cholesky on 512 unknowns, against the accuracy its blocks were compressed atA leaf of 16 on 512 unknowns means 98 truncations inside the factorisation, and the ranks of its blocks run to 4, 7, 10, 13, 16 across the sweep. The residual tracks the tolerance at a slope of 1.031 over ten decades: the knob does what a knob should, and nothing in the approximate arithmetic bends it. The ratio to the representation's own error stays at 0.91, 1.71, 0.86, 0.89, 0.87 — one excursion above one, at 10⁻⁴, where the tolerance sits between two singular values of one block and the rank it takes is a rounding of a decision rather than a measurement.10⁻¹¹10⁻⁹10⁻⁷10⁻⁵10⁻³10⁻¹10⁻¹²10⁻¹⁰10⁻⁸10⁻⁶10⁻⁴10⁻²ε the blocks were compressed at‖A − LLᵀ‖ ⁄ ‖A‖the factorisationthe knob, on a factorisationtruncations98residual at 10⁻⁴2.3·10⁻⁵residual at 10⁻¹⁰9.5·10⁻¹²slope1worst ratio to the representation1.7ten decadesand a slope of one
Fig. 5 The same leaf on 512 unknowns: ninety-eight truncations, ranks 4, 7, 10, 13, 16, residual 0.0017, 2.3·10⁻⁵, 1.4·10⁻⁷, 1.0·10⁻⁹ and 9.5·10⁻¹² at a slope of 1.031. The truncation count has gone from ten to ninety-eight across these two figures and the residual has moved by under a tenth at every tolerance.

The honest reading of the grid is therefore not that two knobs partition every quantity between them. It is that three quantities are in play and only two of them are set by a knob. The leaf sets the shape of the tree and the truncation count; the tolerance sets the accuracy; and the ranks are set by the matrix, with the tolerance choosing where along the matrix’s own decay to cut. A code that treats the rank as something it can configure has misread which of the three it is holding — and that misreading is exactly what the refusal at the foot of this essay is about.

The constant that makes the tolerance quotable

A slope of one is necessary and it is not sufficient. It says a decade in buys a decade out; it says nothing about where the line sits, and a caller asking for a residual of 10⁻⁸ needs the offset as well as the slope.

Measured over all sixty configurations, the residual is smaller than the tolerance at every one of them, and the ratio between them lies between 0.049 and 0.24. That is a factor of five across the whole grid, which is not a calibration. It becomes one when the two ends are separated. At the three tightest tolerances the ratio lies between 0.077 and 0.15, a factor of two; at the two tightest it lies between 0.091 and 0.12, a factor of 1.3 across a twelve-fold spread of trees and a fourfold spread of sizes.

So the rule a code can act on is narrow and it is real: at the accuracies a formatted factorisation is actually run at, ask for a residual, divide by ten, and compress to that. It is good to thirty per cent at 10⁻⁸ and 10⁻¹⁰, and it degrades to a factor of five by the time the tolerance has been loosened to 10⁻⁴ and 10⁻². That degradation is not a defect of the format. A loose tolerance cuts a block after two or three columns, where the decay is steep and the discarded tail is a coarse function of where the cut fell, so the constant is being read off a sequence with almost no terms in it — the same reason a threshold laid across a decaying sequence becomes a decision rather than a measurement when there are few terms to decide between.

This is what the essay’s title is about. A tolerance that is only related to the answer by a slope is a parameter of the method; a tolerance related to the answer by a slope and a constant is a parameter of the problem, quoted in the units the caller asked in. That is the property a stopping rule written in the units of the answer has by construction and that this one has by measurement, and it is the whole difference between a number a code can set and a number a code has to tune.

Where the leaf does reach into the accuracy

The separation above is stated as though it were exact, and there is one place on the grid where it is not.

Each figure carries a second comparison alongside the residual: the error the representation itself commits, before any factorisation is performed. The ratio of the one to the other is the honest measure of what the approximate arithmetic cost, and over most of the grid it sits below one — the factorisation’s residual is smaller than the error of the object it factorised, which is the result the accumulation question turned on.

At a tolerance of 10⁻⁴ that ratio rises above one, and it does so in a pattern the single-size sweep could not have read. At a leaf of 8 it reaches 1.73, 1.78 and 1.67 at the three sizes; at a leaf of 16 it reaches 1.79, 1.83 and 1.71. At a leaf of 32 and at a leaf of 64 it never exceeds 1.00 at any size and at any tolerance on the sweep.

The residual of a formatted Cholesky on 512 unknowns, against the accuracy its blocks were compressed atA leaf of 32 on 512 unknowns means 34 truncations inside the factorisation, and the ranks of its blocks run to 4, 7, 10, 13, 16 across the sweep. The residual tracks the tolerance at a slope of 1.032 over ten decades: the knob does what a knob should, and nothing in the approximate arithmetic bends it. The ratio to the representation's own error stays at 0.92, 0.89, 0.90, 0.89, 0.89 — one excursion above one, at 0.01, where the tolerance sits between two singular values of one block and the rank it takes is a rounding of a decision rather than a measurement.10⁻¹¹10⁻⁹10⁻⁷10⁻⁵10⁻³10⁻¹10⁻¹²10⁻¹⁰10⁻⁸10⁻⁶10⁻⁴10⁻²ε the blocks were compressed at‖A − LLᵀ‖ ⁄ ‖A‖the factorisationthe knob, on a factorisationtruncations34residual at 10⁻⁴1.2·10⁻⁵residual at 10⁻¹⁰9.2·10⁻¹²slope1worst ratio to the representation0.92ten decadesand a slope of one
Fig. 6 Thirty-four truncations at 512 unknowns, with ratios to the representation’s own error of 0.92, 0.89, 0.90, 0.89 and 0.89 across the sweep. The same count of thirty-four at 256 unknowns and a leaf of 16 reaches 1.83 at 10⁻⁴.

That last comparison is the one the grid exists for. Thirty-four truncations at a leaf of 16 produce the excursion and thirty-four truncations at a leaf of 32 do not, so the excursion is not a function of the count. Nor is it a function of the size: it appears at all three sizes at the small leaves and at none of them at the large ones. It belongs to the leaf, which is the one parameter this essay has so far argued has no accuracy term in it.

What it is, is what the figure names: at 10⁻⁴ the tolerance falls between two singular values of one block, and the rank it selects is a rounding of a decision rather than a measurement of a decay. A deeper tree has more blocks, more of them small, and a small block’s spectrum has few terms and coarse gaps — so the chance that the tolerance lands inside a gap rather than on a slope rises with the depth. It appears at a looser tolerance than the ones on either side of it for the same reason, because that is where the gaps are wide enough to fall into.

Two things bound how much this matters. It is a movement in the ratio column and not in the residual column: at 10⁻⁴ and a leaf of 16 the residual is 2.2·10⁻⁵, comfortably inside the tolerance asked for, and the calibration above is untouched. And it is bounded by a factor of two rather than by nothing, across a range of trees that spans the whole grid. The correct statement is therefore not that the leaf is accuracy-neutral but that it is neutral in the residual and not in the margin — a code quoting a residual is safe, and a code claiming that a formatted factorisation is never worse than the representation it was performed on has to say at which leaf, which is a distinction a backward error stated as a bound makes routinely and this field had not needed until now.

The same question one level down

Everything above is about a factorisation. The representation underneath it has the same knob and a different constant, and the two are worth putting side by side because the direction of the offset is not the one a reader would guess.

The accuracy asked for, against the accuracy obtained, for two kernels on one partitionε is applied to each block against that block's own largest singular value; the error is then reported against the whole matrix's norm. The two are not the same number and the dashed diagonal is where they would be. For 1/r the obtained error is 36 times smaller than the tolerance asked for; for log r on the identical partition it is 30 times smaller. Both curves are straight and parallel to the diagonal, so the knob does what a knob should — a decade in buys a decade out — and neither of them sits on it. The factor is how much of the matrix's mass lives off the diagonal, which is a property of the kernel; it is not the number of blocks, which is 46 and would put the curves on the other side of the diagonal.10⁻¹³10⁻¹⁰10⁻⁷10⁻⁴10⁻¹10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³ε asked for, per block‖A − A_H‖ ⁄ ‖A‖ obtainedasked = obtainedlog r1 ⁄ rtwo numbers, not one1/r, obtained at 10⁻⁸2.8·10⁻¹⁰log r, obtained at 10⁻⁸3.3·10⁻¹⁰1/r, obtained ⁄ asked0.028log r, obtained ⁄ asked0.033blocks in the partition46the tolerance is per blockand the error is per matrix
Fig. 7 The tolerance asked for against the error obtained, for two kernels on one partition of 128 unknowns and 46 blocks. For 1/r the obtained error is 0.028 of the tolerance and for log r on the identical partition it is 0.033 — both far below the dashed diagonal.

Compressing a matrix and never factorising it, the obtained error is 0.028 of the tolerance for a 1/r kernel and 0.033 for log r on the identical partition — a factor of about thirty below what was asked. Factorising inside the format, the residual is about a tenth of the tolerance. The factorisation is therefore closer to the tolerance than the representation is, by a factor of three, and it should be: the residual is measured after the approximate arithmetic has been performed, and the representation’s number is measured before it.

The two constants are not the same object and the difference between them is not a defect in either. A per-block tolerance is applied against each block’s own largest singular value and the error is reported against the whole matrix’s norm, so the offset is how much of the matrix’s mass lives off the diagonal, which is a property of the kernel rather than of the partition — which is what makes the two curves in that figure differ at all, since the partition is the same one. What the pair establishes is that both numbers are constants: neither of them drifts across ten decades, so both can be quoted, and neither of them is one, so neither can be omitted. That is the answer to the question a bare tolerance leaves open, asked twice on the same object.

The refusal

The claim under test is one that reads as a definition and is not: that a truncation returns the rank it was asked for.

It is fed a rank-three object and asked to truncate it to rank eight. The result has rank three, because a truncation can only discard, and the assertion that it has rank eight fails. That is the same shape as the mistake this essay is written against. A rank is not a setting; it is what the object turns out to have once a tolerance has been laid across it, and asking for more columns than the decay contains produces nothing — exactly as asking for a different leaf size produced 4, 7, 9, 12, 14 four times over.

The refusal is published rather than asserted in a sentence because it is the boundary between the two readings of the format’s parameters. A code that believes the rank is configurable will size its storage from the number it asked for; a code that has run this will size it from the number the matrix returns, which is the number all twelve trees agree on at a given size and tolerance, and which is the only one that appears in a residual.

What follows for anything that uses the format

A leaf size can be chosen from the machine. It sets the depth, the truncation count and the size of the dense fringe, and it does not set the ranks or the residual. The trade is storage against locality, and it is the same kind of decision as a block size fitted to a fast memory — a machine parameter with no accuracy term, once the margin caveat above is on the record.

A tolerance can be quoted to a caller. Slope 1.022 to 1.046 across twelve trees, and a constant of about a tenth at the accuracies that matter. That is what makes it a knob rather than a tuning parameter, and it is what the site’s rule about printing a residual under every decomposition is for.

A storage estimate needs the size and not the leaf. The ranks moved with the number of unknowns and stood still under the leaf, so the estimate scales with the problem, not with the partition — which is the other half of what an admissibility test buys, where a well-separated block’s rank does not notice the sampling either.

And the crossover argument is unaffected by either. Whether the format is worth using at all is a question of where the count of multiplications crosses a dense factorisation’s, and neither of these two parameters moves the accuracy that argument is conducted at. That is the practical value of a clean separation: the cost argument and the accuracy argument can be settled independently, on the same object, without either of them having to quote the other’s parameter.

What links here

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

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.

AdmissibilityBackward errorCholesky factorisationFormatted arithmeticHierarchical matrixLow-rank approximationOff-diagonal rankRecompressionResidualTruncated SVD