Numerical Linear Algebra
where the algebra stops being the arithmetic
Theme

The thread: Where the algebra stops being the arithmetic — page 7

Essays 145 to 168 of the 262 on this theme, in the same order.
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 Eigenvalues, singular values, rank

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.

 half a stepthe 32 values of one block, in the order they arrivestep 0.0625one scale, thirty-two valuesoctaves inside the block1.9entries rounded to zero0worst error over its bound131 levels either side of zerothe largest entry chose the step The arithmetic underneath

One exponent for thirty-two numbers

Share the exponent across a block and the cost per value drops from eight bits to 6.25, and the accuracy improves — up to about three octaves of spread inside a block. Past that a single outlier deletes the thirty-one values beside it, and the 2-norm barely notices.

 1357911131510⁻¹⁸10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³1vertexrelative error against the exact answerthree routes, one numberworst, by elimination4.8·10⁻¹⁶worst, by spectrum10·10⁻¹⁶solve residual1.5·10⁻¹⁵Foster's sum15n − 115difference3.6·10⁻¹⁵the answer is a ratio of integersso the error is known The matrix that is a graph

A distance computed by a solve

Effective resistance is the one quantity in this field with no combinatorial route to it — it is defined by a linear system. On a small unweighted graph the answer is a ratio of two integers, so for once the error is known rather than estimated, and every resistance in a graph has to add up to a number fixed in advance.

 048121610⁻⁸10⁻⁷10⁻⁶10⁻⁵extra columns in the sample, prelative compression errorthe best rank-8 representationfive draws, mean and rangea band, not a lineexcess at p = 012excess at p = 83spread at p = 00.44spread at p = 80.12the optimum of this rank7.8·10⁻⁸one seed shows the meanand five show the risk Randomised, and the guarantee that changes kind

What a single draw cannot report

With no oversampling the construction's error is 11.6 to 50 times the best representation of the same rank, and the spread across five seeds runs from 16 to 146 per cent of the mean. Sixteen extra columns bring the excess to between 2.0 and 2.9 at every rank measured and the spread to between 4 and 20 per cent — and the second number is the one a single run cannot report and the one that decides whether the first is a measurement.

 A6×6B6×6I ⊗ A + Bᵀ ⊗ I36×36one equation, two objectsentries in A and B72entries in the coefficient matrix1296two routes, relative gap1.7·10⁻¹⁶‖AX + XB − C‖/‖C‖1.4·10⁻¹⁶the small squares are the problemand the large one is the notation Structure, and the solver that cannot see it

The elimination the matrix does not need

The Kronecker form of AX + XB = C is dismissed with a hundred million entries and (2/3)n⁶ operations. Both price an elimination, and after the reduction both routes take, the matrix has exactly n³ nonzeros, none of them above the block diagonal, and nothing left to eliminate.

 234567891010⁻⁴10⁻²1ADI steps krelative error, and max|r|²max|r|²the measured errora bound that is an identitysteps drawn9ratio, worst0.77spread of the ratio1error at k = 103.9·10⁻⁴its residual1.1·10⁻⁴four orders of errorand one constant Reduction, and what a model is for

Bracketing an error nobody can measure

The error of a low-rank Gramian factor is the one quantity a caller cannot compute, because computing it needs the Gramian the factor exists to avoid forming. Two numbers that can be computed sit either side of it — a rational factor known before the run, and a residual known after — and they stay a factor of four apart across a fourfold change of size.

 123456024681012moments Krank of the block Hankel12 eigenvalues insidethe ceiling K·ℓthe ranka ceiling with a knob on itprobes2eigenvalues inside12rank at K = 12rank at K = 612solves, at every K512one ceiling per probeand K of them per moment The eigenvalue problem that is not linear

A ceiling with a knob on it

A contour method returns at most as many eigenvalues as its probe block has columns, and the object that comes back does not distinguish that from having found everything. One line of the derivation multiplies the ceiling by a number the caller chooses, and it costs no extra solves at all.

 point · strong in y0.951point · strong in x0.956y-line · strong in y0.037y-line · strong in x0.967semi-y · strong in y0.107semi-y · strong in x0.9661.00 — no convergenceone problem, seen from two sidesy-line, aimed0.037y-line, turned sideways0.97the ratio2631×31 grid, twelve V-cyclesthe direction is in the code, not in the problem Iterating, instead of factorising

Coarsening in one direction only

Leave the smoother alone and halve only the strong direction, and the smoothing factor is 0.3340 — identical to line relaxation's, at every anisotropy and every weight, to twelve digits. The convergence factors are then a factor of three apart, and at 45° both repairs fail outright.

 -13-11-9-7-5-3-102468log₁₀ of the accuracy asked forlargest train rank0.30 of rank per decadea cost that is typed inentries216stored at 10⁻²90stored at 10⁻¹²288rank per decade0.3worst error ⁄ tolerance0.83the storage is chosena constant of rank a decade When the index is a tuple

The digit that costs more than the tensor

Ask a three-index reciprocal tensor on six points a side for seven digits and its train is 288 numbers against 216 entries. The break-even rank is n − 1 at all four grids measured, and a train that reaches it fits with exactly n numbers to spare.

 unpreconditionedbest step20best error0.14Tikhonov's best0.14cosine, truncatedbest step5best error0.14directions preconditioned25010203040506010⁻¹110¹steprelative errorTikhonov's best: 0.1405plain CGLSpreconditionedthe blur approximated by a cosine transformtruncation leaves it alone Methods that were designed apart

What a cheap preconditioner has to leave alone

A blur approximated by a matrix the cosine transform diagonalises agrees with the operator everywhere but its first and last seven rows. Made invertible by a shift, as the exact preconditioner was, it never reaches the unpreconditioned run's floor — at α = 10⁻³ its best iterate is 0.749 against 0.143. Made invertible by leaving every eigenvalue below τ alone, it reaches 0.141 in five steps instead of twenty, and the smallest τ that keeps the floor sits at a third to a half of the Tikhonov oracle's λ at three noise levels.

 122028364452606876849210010⁻¹110¹grid points nrelative error against the continuous signalbest with no λ: n = 26no λeach grid's best λTikhonov, 64 points0.1% noise, 16 drawsbest grid with no λ26its error0.1396 points with its λ0.12best λ on fine grids0.0056each point is a median over the same drawsλ belongs to the problem, not to the grid Regularisation, and the answer that is chosen

Where the grid hands over to λ

An unregularised solve on a coarse grid comes within a tenth of the best Tikhonov answer on a fine one, and the pair of a grid and a λ was left unmeasured. Measured, the two do not trade. On every grid up to the best unregularised one no λ helps at all. On every grid of 40 points and more the best λ is the same to within a quarter of a decade — 3.2·10⁻² at 1% noise per sample, 10⁻³ at 0.01% — and the 96-point grid with it beats the best coarse grid by 7, 9 and 13 per cent. The grids between the two, given their own λ, land between them.

 least possible76minimum degree76worst of 300 tie-breaks81reverse Cuthill–McKee85natural order9920 unknownstie-breaks at the least0.44minimum degree ÷ least1every order searched, by the set already eliminatedthe greedy rule found the minimum Sparsity, and what elimination costs

The least fill there is

Finding the elimination order with the least fill is NP-hard, and that is a statement about the hardest graph and the largest size. On a graph of twenty vertices every one of the 20! orders can be searched at once, through the million sets of vertices already eliminated, and the least fill is a number. On the 4×4, 4×5 and 3×7 grids minimum degree finds it exactly. On eighty random sparse graphs of eighteen vertices it finds it on 53 and misses by at most 7.6 per cent, and on every one of the eighty some breaking of its ties finds it.

 0246810121416110¹10²10³10⁴10⁵10⁶length of the interval Tgrowth factorpartialrookcompleteh = 0.3, 102 unknownspartial, T = 15.01.3·10⁵e^(5T/6)/21.3·10⁵complete2largest κ8.3a boundary-value problem, not a constructionκ single-digit throughout Elimination, and the swap

The growth a boundary-value problem supplies

Large growth under partial pivoting is usually said to need a matrix built for it. A two-point boundary-value problem solved by multiple shooting supplies one without being asked: its growth factor is e^(5T/6)/2 to four figures — 1.1·10⁴ at an interval of 12, 1.3·10⁵ at 15 — on a matrix whose condition number never exceeds 8.3, while rook and complete pivoting keep it below 2. And it is the finer shooting grid that grows: below a step of 0.3397 the choice partial pivoting makes turns on one entry against one, and above it there is no growth at all.

 10⁴1ceiling on the largest intermediate, numbers heldarithmetic ÷ the unconstrained cheapest orderno evaluation order fits below thisminimising the peak insteadthe unconstrained cheapest ordera ceiling, not a targetfree order, peak9504leanest possible2106cheapest that fits1leanest order costs1.6the difference1.6largest input2106the same memoryat two prices Where the flop count stopped predicting the time

A ceiling is not a target

Asked to hold the smallest possible intermediate, an evaluation order costs a median of 1.55 times the cheapest order's arithmetic. Asked to hold no more than that same amount, it costs 1.15. The two answers hold exactly the same number of numbers, and on one network they are 5.17 times apart in work.

 0204060801001share of shufflings, per centarithmetic ÷ the sorted order'sthe sorted order, and the nearest-neighbour paththe order is a free variabledrift between neighbours0.04median penalty1.5worst penalty1.8factorisations, sorted4worst shuffle's16greedy ÷ sorted1the same sixteen problemsat two prices When the problem arrives again

The order a batch arrives in

Sixteen problems over a parameter, solved in the order the loop produced them, cost a median of 1.54 times what the same sixteen cost sorted, and 2.80 times at the worst shuffling. A nearest-neighbour path computed from the parameter values alone recovers the sorted cost exactly, at every drift and every shuffle.

 13579111310⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²10¹index ksingular valueσ₁ · 10⁻¹⁴, the thresholdone matrix, three rankspartitionings7lowest rank10highest rank12threshold7.3·10⁻¹⁴σ₁7.3the curves separate in the noiseand the threshold is drawn through it The answer that depends on the machine

A rank that depends on the thread count

One 60 × 14 matrix, one threshold, seven partitionings of the inner products that build its Gram matrix — and numerical ranks of 12, 12, 12, 10, 10, 11 and 11. Not a digit of an answer: the number of columns a model built from this matrix would have.

 19172533414911.251.51.752indexsample ÷ originalkept, and not keptedges kept344of1225eigenvalue ratio, low0.43high1.7degree ratio, low0.56high1.5diameter, after3filled dots are eigenvaluesopen dots are degrees The matrix that is a graph

A graph with a tenth of the edges

Keeping 344 of 1,225 edges, sampled by effective resistance and reweighted, preserves every eigenvalue of the Laplacian to within a factor of 1.7. It preserves no degree — half of them are wrong by more than a third — and it takes the diameter from one to three.

 0246810121410¹⁵10¹⁹10²³10²⁷10³¹10³⁵10³⁹10⁴³k|CAᵏB|6 exact zerosa definition that will not computestates16exact zeros6‖A‖₂, measured1146last step-to-step ratio1330range across the run2·10²⁹zero for the travel timethen ‖A‖₂ᵏ Reduction, and what a model is for

The definition asks for more of what defeats it

The rank of a p × p Hankel matrix of Markov parameters resolves a degree of p, and it needs 2p parameters to do it. Those parameters grow like the norm of the state matrix raised to their index, so the count that buys resolution is the same count that buys dynamic range. One model, four run lengths, and a spread that runs from 10¹¹ to 10¹⁰².

 entries rounded to zero, of 320block, as given310block, sorted by size22E4M3, either order0median entry's relative errorblock, as given1block, sorted0.00321E4M30.022the same numbers, three waysdeleted, as given310deleted, sorted22deleted, per-element08-bit significands, blocks of 32sorting is free and changes no value The arithmetic underneath

A bit buys an octave

The outlier a block survives is exactly two raised to its significand width — 8 at three bits, 32 at five, 128 at seven, 512 at nine. Each extra bit doubles the range the block tolerates and halves the ordinary entry's error. Reordering the same numbers buys every octave at once and costs nothing.

 123456110²10⁴10⁶moments Kκ of the block Hankelintegrated in zin (z − c)/ρone division per pointradius6unscaled, K = 62.3·10⁶scaled, K = 61750the factor between1291the ceiling rises by Kand so did the conditioning The eigenvalue problem that is not linear

The conditioning that rises with the ceiling

Higher moments multiply a contour method's ceiling by K and grade its block Hankel over ρ to the 2K, so the two knobs are the same knob. One division per quadrature point separates them, and the measurement of what it is worth grows from twenty to twenty thousand.

 -12-10-8-6-4-2010⁻¹⁷10⁻¹³10⁻⁹10⁻⁵10⁻¹10³10⁷10¹¹10¹⁵log₁₀ μcondition number, and relative errorκ₂, augmentedcomponentwise, condensedcomponentwise, augmentederror, augmentedone matrix, two numbersκ₂ at μ = 10⁻¹²3·10¹³componentwise, same matrix13their ratio2.3·10¹²measured error9.4·10⁻¹⁶every library prints the top lineand the error obeys the third Two errors, and whose fault they are

Two condition numbers of one matrix

κ₂ is a worst case over perturbations of a given norm, and a normwise perturbation may put its whole budget on the smallest entry. The componentwise number is a worst case over perturbations proportional to the entries, which is what a backward-stable factorisation actually makes. On one matrix they are 3·10¹³ and 13.3, and the error obeys the second.

 10¹10²110¹rank of the trains at the moment of evaluationcost ÷ the cost of this rank's own best ordereach rank's own best orderthe order compiled at rank 4, carriedthe smallest-result rule, recomputedthe cheapest-product rule, recomputedone plan, eight rankscompiled at rank4worst carried8carried at rank 2568fresh rule at 2562.9carried at rank 21fresh rule at 21.1the search was exhaustiveand the dimensions moved Where the flop count stopped predicting the time

The plan that was right at rank four

An evaluation order is chosen once and paid for thousands of times, and the dimensions it was chosen at are not the dimensions it runs at. Compiled at rank four and run at rank 256 it costs 8.01 times the order that rank deserves; compiled at 256 and run at 2 it costs 301 times. A one-line rule recomputed on arrival costs 2.92 and 1.11.

 38434853586368110¹10²10³10⁴e, where the skew is 2ᵉ + 12345orthogonality defect of the basis returneda double holds k exactlyan orthogonal basisthe invariant sees none of itfloat, at 2^521float, at 2^531.4float, at 2^701.2·10⁴residue at 2^701.2·10⁴exact, everywhere1determinant, always1every step was unimodularand the answer is 12,345 times worse Exact arithmetic, and what it costs instead

The knob and the rounding

The Lovász parameter is the number a lattice reduction is specified by, and moving it from 0.50 to 0.99 strengthens the proved bound from 16.00 to 1.83, costs 91 per cent more steps, and returns a basis with the same orthogonality defect. The one rounding nobody writes down decides everything: above 2⁵³ the reduction returns a basis 12,345 times worse than it should, with the determinant invariant equal to one throughout.

 01234567110¹10²levelstored entries per rowalgebraicgeometrictwo complexities, one hierarchygrid complexity3.1operator complexity18geometric, for comparison1.5a ring plus random chords — nothing is anywherelevel 6 is 100% dense Iterating, instead of factorising

A hierarchy with no grid behind it

On a graph Laplacian the algebraic V-cycle converges at 0.199 a cycle, its grid complexity is an unremarkable 3.05, and its operator complexity is 17.7 — one level of forty-one unknowns is entirely dense. The number people quote is the one that does not measure the work.

All themes