Numerical Linear Algebra
where the algebra stops being the arithmetic
Theme

The thread: Two routes to a number — page 6

Essays 121 to 144 of the 228 on this theme, in the same order.
-18-15-12-9-6-300log₁₀ ‖PKPᵀ − LDLᵀ‖ / ‖K‖share of orderingsbestworstexistence and stabilityfactorise1unregularised0.67worst growth5·10⁷growth × γ0.5the ordering is free to chooseand not free of consequence The matrix a constraint makes

The perturbation that does the work

A saddle-point matrix made quasi-definite is perturbed in both blocks, and the laws measured for it moved both together. Moved apart, the laws all belong to one block. The zero block's perturbation γ decides whether every ordering factorises, sets the worst ordering's growth at 0.51/γ, and costs the answer 1,451 per unit — the reciprocal of the smallest eigenvalue of AH⁻¹Aᵀ to three figures. The perturbation of H moves none of the first two and costs 19 per unit. Refinement removes each block's perturbation at the rate its own Schur complement sets, so γ's limit sits fifty times nearer than δ's.

 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.

 the grid the operator came froman edge is a coupling the matrix calls strongkeptinterpolatedwhat the entries decidedcoupling ratio, x against y0.001strong couplings across x0strong couplings along y110rows kept or dropped whole1111×11 grid, θ = 0.25the coarse grid, from the matrix alone Iterating, instead of factorising

The coarse grid the matrix chooses

Given a tridiagonal matrix and no information about a grid, the coarsening keeps every other point and derives the weights ½, 1, ½ — the operators the geometric method was handed. Given the anisotropic operator, it discovers semi-coarsening, in the right direction, without a coordinate.

 10⁻³10⁻²10⁻¹110⁻¹⁹10⁻¹⁷10⁻¹⁵10⁻¹³10⁻¹¹10⁻⁹1 − ρbackward errorLevinson, Toeplitzelimination, ToeplitzLevinson, any kindelimination, any kindat ρ = 0.999Levinson, any perturbation2.9·10⁻¹⁷Levinson, Toeplitz only6·10⁻¹²the ratio between them2·10⁵diagonal defect of the first0.97the number that is reportedand the number that is asked about Structure, and the solver that cannot see it

The number that cannot rank them

Levinson and Gaussian elimination are indistinguishable on the backward error a library reports — every one of ninety-six measurements between 1.16·10⁻¹⁷ and 5.73·10⁻¹⁷. The structured backward error separates them by up to a hundredfold, in whichever direction the point happens to give. Only the forward error ranks them, and only because this family's exact answer is known.

 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¹⁰².

 1234567891010⁻⁵⁵10⁻⁴⁹10⁻⁴³10⁻³⁷10⁻³¹10⁻²⁵10⁻¹⁹10⁻¹³10⁻⁷10⁻¹cyclesizethe residualrepairedthe reported boundwhat the stopping rule readsreported at the last cycle9.4·10⁻⁴¹the residual there5.7·10⁻⁵with the border recomputed3.8·10⁻¹⁴products, cheap and repaired9a bound with nothing under itand one product a cycle to fix it Eigenvalues, singular values, rank

Keeping the vectors, and losing the bound

Thick restarting keeps the Ritz vectors instead of filtering the starting vector — the same eigenvalues for a third of the products with A. Its residual bound reaches 9.4·10⁻⁴¹ while the residual it bounds sits at 5.7·10⁻⁵, and the eigenvalues are correct to 4.3·10⁻¹⁴ the whole time, so nothing reports it.

 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.

 8×8, κ = 10⁸direction at ε = 10⁻², ‖QᵀQ − I‖1.5·10⁻¹⁵scalar at ε = 10⁻², ‖QᵀQ − I‖0.065direction at ε = 10⁻², residual0.00510⁻¹⁴10⁻¹²10⁻¹⁰10⁻⁸10⁻⁶10⁻⁴10⁻²10⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²relative perturbation εdeparture‖QᵀQ − I‖, scalar‖A − QR‖/‖A‖, scalar‖A − QR‖/‖A‖, direction‖QᵀQ − I‖, directionflat line: a different reflection is still a reflectionsloped line: a non-reflection is not Orthogonality, measured

One number that has to be right

Householder's orthogonality was called structural: a reflection is built from a unit vector, so rounding the vector names a different reflection rather than a broken one. Tested by breaking it, the claim is narrower and sharper. Perturb every component of the reflector by a relative 10⁻², and ‖QᵀQ − I‖ stays at 1.5·10⁻¹⁵ while the factorisation moves to 5·10⁻³. Perturb the one stored scalar by the same amount and ‖QᵀQ − I‖ is 6.5·10⁻². The structure is one degree of freedom, and the departure is four times its relative error.

 all 2^16 subsets searchedboth minima, one order0.97minimum degree optimal0.82worst fill overshoot1.1worst work overshoot1.111.021.041.0611.041.081.121.16fill ÷ the least filloperations ÷ the leastovershoot squaredon the diagonal: the two overshoots are equalon the upper curve: the arithmetic overshoot is the square Sparsity, and what elimination costs

Two minima that are one minimum

The order that decides the memory found the operation count behaving like the square of the fill, which leaves room for an order with slightly more fill but a shorter heaviest column to do less arithmetic. Searched exactly over every elimination order on forty graphs, that order does not exist: one order attains both minima on thirty-nine of forty, and on the fortieth the least-fill order's arithmetic is 1.0099 times the least. Minimum degree attains both on the same thirty-three graphs and neither on the same seven.

 12345678910111210⁻¹⁵10⁻¹³10⁻¹¹10⁻⁹10⁻⁷10⁻⁵10⁻³10⁻¹degree of the approximantworst relative error on the target setchosen by the residualdegree reached12error there10⁻¹²even support, same degree1.1·10⁻⁷advantage702nearest support point0.096clustering ratio1.5filled: points chosen by the erroropen: points spread evenly The eigenvalue problem that is not linear

The points the algorithm chose

A rational approximant whose support points are picked by its own residual clusters geometrically at a branch point nobody named — recovering by measurement the rule a hand-built approximant is given. At degree ten it is seven hundred and fifty times more accurate than the same form with its points spread evenly.

 ‖H − Hₙ₋₁‖∞ ÷ σₙpositive residues2.0000mixed residues2.0000clustered poles2.0000four poles2.0000wide spectrum2.0000‖ΔH‖ 3.89·10⁻⁴ σ 1.94·10⁻⁴‖ΔH‖ 3.92·10⁻⁴ σ 1.96·10⁻⁴‖ΔH‖ 0.00502 σ 0.00251‖ΔH‖ 8.41·10⁻⁴ σ 4.21·10⁻⁴‖ΔH‖ 4.6·10⁻⁴ σ 2.3·10⁻⁴2four systems, four ratios, one number: 2.0000one state, one equalityand the sum is that, repeated Reduction, and what a model is for

The state that is removed is not a mode

Balanced truncation removes one state and pays exactly twice one Hankel singular value. The natural reading is that the state removed is the model's least important mode and that the σ is that mode's own size — and on five systems that reading over-estimates by between 1.21 and 7.46, never once under.

 12345678110¹partial orders kept at each pairingexcess over the exhaustive orderthe exhaustive orderdashes above: each rule's worst casepair by cheapest productpair by smallest resultone knob, two rulesby cost, width 11.5by cost, width 21.1by size, width 11both, width 81by size, worse wider24by cost, worse wider3the width repairs the bad ruleand does not improve the good one Where the flop count stopped predicting the time

The search that got worse as it widened

Keeping two candidate orders instead of one removes a third of the cheapest-product rule's excess, and by eight it has removed all of it — the two greedy rules this field separated become the same rule. On twenty-four of sixty networks a wider search returns a dearer order than a narrower one, and it is the better rule that it more often makes worse.

 unpreconditioneddimension a step1.2best iterate carries24at step20τ = 10⁻³, past the edgedimension a step3.6best iterate carries33its error0.4837111519232731353910⁻¹10⁻⁰.⁵1effective dimensionrelative errorthe arc: dimension 8.0 to 24.0none10⁻⁰.⁵10⁻¹10⁻¹.⁵10⁻²10⁻².⁵10⁻³large dots: each run's own best iteratethe answer sits at one effective dimension Methods that were designed apart

The parameter neither knob is

A preconditioned run has a cutoff and a step count, and neither is the regularisation parameter. The parameter is the effective dimension of the iterate: every cutoff that works puts its own best at 23.7 to 24.3 of it, where the unpreconditioned run's best sits at 23.2, and what the cutoff buys is the rate — 1.27 of it a step with no preconditioner and 3.53 with one. The edge is where a single stride is longer than the distance left.

 what each reachesbest pair0.0034‖x‖ alone0.0048‖L₁x‖ alone0.0035λ₁, on ‖x‖λ₂, on ‖L₁x‖ — log₁₀0-6-5-4-3-2-100-6-5-4-3-2-10outlined row and column: one penalty switched offred cell: the best pair, 0.00343darker is worse, on a logarithmic scalethe optimum sits on or beside an edge Regularisation, and the answer that is chosen

A second penalty is not a second parameter

Penalise ‖x‖ and ‖L₁x‖ at once and there are two λ to choose. Over ten draws on five signals the best pair beats the better single penalty by between 0.00% and 3.1%, and one of its two parameters is exactly zero on 30 to 70 per cent of draws. Choosing the wrong one of the two costs up to 54%. The surface is a choice between two curves with a knob nobody needs.

 10×10, ±6determinant, bits27Hadamard's bound36median slack10024680510152025303540elimination stepbitsHadamard, on the minorsthe intermediatesevery intermediate is a minor of the originalso a bound on the minors bounds the run Exact arithmetic, and what it costs instead

A bound on every intermediate at once

Fraction-free elimination's intermediates are minors of the original, which is a theorem about exactness. It is also a bound: Hadamard's inequality applies to every minor, so one inequality bounds the whole run before it starts. The bound on the k-th step is the one on (k+1)×(k+1) minors, not the one on the whole matrix — and on a 10×10 with entries in ±6 the difference is ten bits, with the run reaching 2.7 bits a step against the bound's 3.4.

 six tolerancesblocks, every tolerance250slope per unit of rank5.1worst residual3.6·10⁻¹⁵234567010203040mean rank of the low-rank blocksnumbers per doubling2 digits4 digits6 digits8 digits10 digits12 digitsthe partition is the same at every toleranceand the rank rises one per two decades Neither sparse nor dense

The partition that does not move

The storage curve's growth rate more than triples between two digits and twelve, and the earlier measurement said so without saying what moves it. There are two candidates and the measurement separates them decisively: the partition is identical at every tolerance — 250 blocks, 156 of them low-rank, 94 dense — and the mean rank rises by exactly one per two decades of accuracy. The slope is 5.11 numbers per unknown per doubling per unit of rank, with a worst residual of 0.03 over six tolerances.

 081624324000.250.50.7511.25angle of the strong direction (degrees)ratio, and convergence factorenergy ratio of oneenergy ratioclassicalminimiserwith the constraintwhere the formula is optimalenergy ratio at 0°1energy ratio at 45°1.2classical rate at 45°0.35the minimiser's rate there0.28a ratio of exactly one while the assumption holdsand a diagnostic when it stops Iterating, instead of factorising

The formula that was already optimal

Ask for the interpolation that minimises the energy of its own columns and the answer is the classical AMG formula — to zero at every row of the one-dimensional Laplacian, and to four digits in two dimensions. On the operator rotated to 45° the two part company, and the gap between them is a diagnostic that needs no reference solution.

 10⁻¹110¹10²10³1234eigenvalue of the preconditioned Laplacian, and the stretchbreadth-firstκ = 34.62randomκ = 55.17random, secondκ = 37.11the bar is the spectrumthe tick is a combinatorial bound The matrix that is a graph

A preconditioner that is a tree

Every eigenvalue of a tree-preconditioned Laplacian is at least one and at most the total stretch — a combinatorial integer with no arithmetic in it. Measured, the bound is two to four times loose, and on a grid the preconditioner makes the conditioning worse by a factor of 1.85 at every size.

 10¹10³10⁵10⁷10⁹10¹¹10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³1κ(A), the conditioning of the fitrelative error in the answerthe saddle-point route, in floating pointdashes: forming AᵀA, then solving exactlythe method of weighting, τ = 10⁸the null-space routethe reference was a methodκ(A)10¹¹κ of what is left4.9·10⁶null space1.9·10⁻⁹saddle point4.6·10⁻⁴forming AᵀA alone3.4·10⁻⁴weighting3·10⁻⁹the damage is in the formingand not in the solving Least squares, and the road not to take

The reference was a method

The optimality conditions of a constrained fit contain AᵀA, so solving them is the road that squares the problem wearing a block structure. At κ(A) = 10¹¹ the route that never forms a cross-product returns 1.89·10⁻⁹ and the route that does returns 4.64·10⁻⁴ — and forming AᵀA and then solving it in exact rationals returns 3.45·10⁻⁴, so nearly all of the loss happens before any elimination begins.

 -7-6-5-4-3-210⁻¹⁶10⁻¹⁴10⁻¹²10⁻¹⁰10⁻⁸10⁻⁶10⁻⁴10⁻²110²10⁴log₁₀ of the twist's distance from landing, in turnsrelative forward errorsimple zero, landed on at 0simple zero, landed on at 0.3double zero, d = 2n = 64, twist 10⁻⁵ of a turn from landingpair: error, κ(A) = 4.1·10⁶5.4·10⁻¹²single: error, κ(A) = 3.8·10⁶3·10⁻⁸double: error, κ(A) = 4.2·10¹²7515the same condition number of the wrapa different number of modes near zero Structure, and the solver that cannot see it

Two near-zeros cost less than one

Solve a well-conditioned tridiagonal matrix through a nearly singular wrap and the correction's accuracy is not set by how singular the wrap is. At κ = 4·10⁷ one wrap returns the answer to 8.8·10⁻¹¹ — better than κ·u — and another, at κ = 3.8·10⁷, returns it to 3.3·10⁻⁶. The difference is how many of its samples sit near the symbol's zeros. A real wrap lands on a conjugate pair, a rank-two correction absorbs the pair exactly, and its two-by-two system has condition number 1.00 — which mattered because that system was solved by Cramer's rule; solved with pivoting, the single landing costs what the pair costs.

 11.21.41.6partial orders kept at each pairingcost over the exhaustive order's cost123468cost so far+ largest group+ closing pairings7 tensors, 60 networkscost so far: networks worse wider24+ largest group: networks worse wider12+ closing pairings: networks worse wider23solid: median · dashed: ninetieth percentilethe width axis is the same knob in every curve Where the flop count stopped predicting the time

A beam ranked on what remains

A beam over contraction orders ranked on cost so far returns a dearer order when it is widened on 24 of 60 networks. Rank it on cost so far plus the largest group still to be paired — a lower bound on what remains, and free — and that falls to 12, eleven of them from the original 24. Rank it on a stronger estimate that is not a bound and the count stays at 23, but only 10 are the same networks: the anomaly has moved, not gone. At nine tensors, where the unranked beam's median at width 8 was worse than at width 1, both rankings make widening pay again.

 10⁻¹101020304050z, the effective noise over thicknesstrials mirrored, %5 points, equal20 points, equal80 points, equal● five precise, 1/σ²■ five weighted ×100per cent mirrored at z = 0.45, 1,000 trialstwenty points1.8eighty points1.3five points9.5five precise, weighted 1/σ²7.9equal noise, five weighted ×1006.8both weighted schemes have five effective pointstwenty and eighty points share one curve Orthogonality, measured

Five precise points are five points

Weighting each sighting by its reliability is the standard form of an attitude or registration fit, and it changes how often the nearest orthogonal matrix comes back as a mirror. Measured, the rate is a function of two numbers: the weighted noise over thickness, and the effective count (Σw)²/Σw². Five points with a tenth of the noise, weighted by 1/σ², carry the information of 515 equal points and mirror like five — 7.9 per cent at a noise ratio where twenty points mirror 1.8 and five mirror 9.5. The √m the earlier measurement left unchecked is right, and it counts what carries the thin direction.

 57911131517192110⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²10¹vertices of the complete graphrelative sizethe algorithm is not the problemsizes counted18worst relative error9.4·10⁻¹⁵first wrong at n =16digits in the answer17spacing of doubles there8largest integer a double holds9·10¹⁵the elimination stays at the rounding levelthe format runs out The matrix that is a graph

A count that comes out of a determinant

The number of spanning trees of a graph is the determinant of its grounded Laplacian, so it is a whole number known in advance. The elimination that computes it is backward stable at every size — and from sixteen vertices the answer is wrong, because the count has seventeen digits and a binary64 has sixteen.

All themes