Numerical Linear Algebra
where the algebra stops being the arithmetic
Concept

Hierarchical matrix — where it appears

A matrix stored as a partition of blocks, most of them held as two thin factors and the rest densely. Its storage is a function of the accuracy asked for rather than of the matrix, which is the property no other format here has.

Named by 11 essays across 5 fields — each of them below, with the objects they name alongside it.

45678903691215log₂ of the points a sidecolumns above 10⁻⁸two intervals that touchtwo that do notthe measurement that is a null resultadmissible, n = 325admissible, n = 2565touching, n = 329touching, n = 25613stored ⁄ dense at largest0.039the rank belongs to the geometryand not to the sampling

The size the rank does not notice

Sample a kernel block at 32, 64, 128 and 256 points a side and it needs five columns, five, five and five. Sample the touching block next to it at the same four sizes and it needs nine, eleven, twelve and thirteen. Same kernel, same accuracy, one number and a logarithm.

hierarchy · off diagonal rank -1012301020304050log₂ of the segment lengthcolumns above 10⁻⁸cos(κr) ⁄ r at κ = 401 ⁄ r, same pointsthe case with no answerq, held fixed0.51/r at L = 0.561/r at L = 86cos(κr)/r at L = 0.512cos(κr)/r at L = 853the geometry did not moveand the rank did

The kernel with nothing to compress

Hold the geometry fixed at q = ½, fix the wavelength, and scale the picture up by sixteen. A smooth kernel needs six columns at every scale. An oscillatory one needs twelve, sixteen, twenty-two, thirty-three, fifty-three, and there is no scale at which it stops.

hierarchy · off diagonal rank 545545545455554545545545545454555555454545545545545455554545545545rows and columns, in the order the points arrivedark: kept dense · light: two thin factors, rank printedthe strong partitionblocks112kept dense46largest rank5numbers stored2.7·10⁴‖A − A_H‖ ⁄ ‖A‖3.4·10⁻¹⁰the picture is decidedbefore a number is read

Which pairs are allowed to be small

A hierarchical representation is a partition of the matrix into blocks, and the rule that produces it reads four numbers per pair of index clusters and not one entry of the matrix. On a 256-square it yields 112 blocks, 66 of them stored as two thin factors, none of rank above five.

hierarchy · admissibility 567891003691215log₂ nlargest rank in the partitionweak: touching pairs compressedstrong: touching pairs refusedthe test bounds a rank and costs storagestrong, blocks250weak, blocks94strong, numbers6.8·10⁴weak, numbers6.1·10⁴strong ⁄ weak1.1the better partitionis the more expensive one

The test that costs what it saves

The partition that refuses to compress a touching pair keeps every rank at five while the other lets them climb from nine to thirteen. It also stores more numbers at every size measured — 67,968 against 61,440 at n = 512 — and which of those two facts matters is a question about how large the problem is going to get.

hierarchy · admissibility numbers stored, 256 × 256clustered, strong27,008clustered, weak24,064the dense matrix65,536shuffled, strong65,536shuffled, weak118,208the same matrix, twiceκ, clustered24κ, shuffled24‖A‖_F, clustered6140‖A‖_F, shuffled6140shuffled weak ⁄ dense1.8the compressibility is in the numberingand the numbering is not in the matrix

The same matrix, numbered twice

One symmetric permutation. The condition number is 24.3948 either way to eight digits and the Frobenius norm is 6.13996414·10³ either way to twelve. The partition that stored 27,008 numbers now finds no admissible pair anywhere and stores all 65,536, and the format that compresses regardless stores 118,208.

hierarchy · admissibility 567891010⁴10⁵10⁶10⁷10⁸10⁹log₂ nmultiplicationsdense factorisation, n³⁄3the recursion, countedwhere the format starts payingratio at n = 641.5ratio at n = 5120.16exponent, first doubling2.1exponent, last doubling1.7backward error1.4·10⁻¹⁰cheaper is a sizenot a property

Where the format starts paying

A hierarchical solve costs 1.48 times a dense factorisation at 64 unknowns and 0.16 times it at 512. The crossover is between 64 and 128, it walks right when the accuracy is tightened, and the exponent between consecutive sizes is 2.13, 1.93, 1.74 — falling towards one and never arriving.

cost · hierarchical solve 567891010¹10²10³10⁴10⁵10⁶log₂ nnumbers the construction touchedentries, n²products with the operatoran operator, applied a few hundred timesproducts at n = 512256entries at n = 5122.6·10⁵per doubling48‖A − A_H‖ ⁄ ‖A‖4·10⁻⁷excess over the compression7.3no entry of the matrixwas ever read

Built from products alone

A 512-square hierarchical representation, at a relative error of 4·10⁻⁷, from 256 applications of an operator that is never assembled. The compression route reads 262,144 entries; this one reads none, and pays for it with a factor of seven against the representation the entries would have given.

randomised · sketching -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

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.

hierarchy · hierarchical solve 0481216202402468101214unknowns on the separatorcolumns above 10⁻⁸the same matrix, renumberedin the separator's own orderthe ordering the geometry hands overseparator 73separator 236renumbered, largest11the block, largest11share of the square stored0.52the fill is totaland it is not independent

The fill that is not independent

Eliminate both halves of a grid and what is left on the separator is 100 per cent nonzero — the sparsity field's result, unchanged. Its off-diagonal block is 11 by 12 and six columns describe it to eight digits. Renumber the separator and the same block needs all eleven.

sparsity · fill 02040608010010⁻¹¹10⁻¹⁰10⁻⁹10⁻⁸truncations performed by the factorisationrelative errorthe representation's own error‖A − LLᵀ‖ ⁄ ‖A‖no decomposition without its residualresidual, 0 truncations2.1·10⁻¹⁰residual, 98 truncations1.1·10⁻⁹representation, deepest1.4·10⁻⁹residual ⁄ representation0.81levels, deepest5a hundred approximate stepsand an exact-looking factorisation

The rounding that was not the problem

A rank-k block plus a rank-k block is a rank-2k block, exactly, so every arithmetic in this format truncates after every addition. A Cholesky performed inside it does ninety-eight of those and its residual is 1.14·10⁻⁹ against a representation error of 1.40·10⁻⁹ — the roundings cost nothing measurable.

hierarchy · recompression 10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³‖A − A_H‖ ⁄ ‖A‖, the representation‖b − Ax‖ ⁄ ‖A‖‖x‖, the solveequala backward error, chosenslope0.99representation at 10⁻⁸1.2·10⁻⁹backward error there1.3·10⁻¹⁰representation at 10⁻¹²1.2·10⁻¹³backward error there2.5·10⁻¹⁴the compression is not an approximationit is a perturbation of the problem

An accuracy that is a backward error

Every backward error on this site is something an algorithm produced and somebody then measured. This one is a line in the program. Solving with a compressed matrix gives a residual that is the compression's own error, at a slope of 1.000 over ten decades, so the knob that sets the storage sets the backward error directly.

error · backward error

Named alongside it

The objects these essays reach for when they reach for this one.

AdmissibilityOff diagonal rankAsymptotic analysisKernel matrixNumerical rankBackward errorBlock methodsCluster treeCondition numberFlop countWoodbury identityElimination order

All concepts