Numerical Linear Algebra
where the algebra stops being the arithmetic
Theme

The thread: Measured, not asserted — page 14

Essays 313 to 327 of the 327 on this theme, in the same order.
33.310⁻⁴10⁻³10⁻²10⁻¹1log₁₀ numbers storeddistance to the dominant eigenvalueArnoldi, linearisedprojected quadraticper number heldstorage, linearised2080storage, second-order1040Ritz values, linearised26Ritz values, second-order52half the storageand twice the approximations Iterating, instead of factorising

A Krylov space for a problem that is not linear

A quadratic eigenvalue problem has no matrix to build a Krylov space out of. The recurrence that builds one anyway stores half as many numbers, returns twice as many Ritz values — and stops being a basis at twenty vectors while the answer it gives keeps improving.

 0102030405010⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²10¹10⁴10⁷decades of gradingworst relative errorthe answer is gonevia BᵀBone-sided Jacobizero-shift QRone axis, four routesBᵀB at the narrowest grading3.6·10⁻¹⁴and at the widest1.9·10⁷Jacobi, worst over the sweep1.5·10⁻¹⁵zero shift, worst1.1·10⁻¹⁵the definition is not a methodand squaring buries what it squares Eigenvalues, singular values, rank

A threshold the matrix does not set

Two numbers come out of a relative-accuracy comparison and they belong to different things. The size of the matrix moves the constant of the routes that never fail, by a factor of 2.7 between n = 4 and n = 10; it does not move the point where the route through BᵀB stops returning an answer, which sits between ten and eleven decades of grading at every size drawn.

 00.250.50.75100.250.50.751xuexact at ε = 0.005exact at ε(1 + Pe)the upwind answeran exact answer to a different question‖upwind(ε) − central(ε(1+Pe))‖/‖·‖0distance to the problem it solves0.026distance to the problem posed0.36the added diffusion is h/2 = 0.001953, whatever ε isso refining removes it Iterating, instead of factorising

A different equation on every grid

Upwinding is the exact discretisation of a convection–diffusion problem with diffusion ε + h/2, entry for entry, at a relative difference of between 0 and 1.26·10⁻¹⁶ on every mesh from 15 points to 511. The equation it is exact for is chosen by the mesh and not by ε — the added diffusion is 0.01563 on a 31-point grid whether ε is 0.2 or 0.001.

 01020304050607010⁻⁶10⁻³110³10⁶10⁹k‖Aᵏ/k!‖‖e^A‖ = 2.6largest term 1.4·10⁷what the series throws awaylargest term1.4·10⁷‖e^A‖2.6digits cancelled away5.4·10⁶error after the sum5.2·10⁻⁹every term is computed correctlyand the sum has lost seven digits Eigenvalues, singular values, rank

The error the method already knows

Summing the exponential's Taylor series throws away a known number of digits, and the number is on the machine while the sum is being formed. The largest term divided by the answer, times the unit roundoff, tracks the relative error that comes out — to within a factor of nine, across fourteen orders of magnitude of it — and nothing reports it.

 10²10⁻⁴10⁻³10⁻²10⁻¹grid points nworst nodal errorupwindtunedcentralthe same tuning, another problemtuned ÷ central at n = 3149tuned ÷ central at n = 12796central's error at the finest grid10·10⁻⁵exact on the problem it was derived fromand harmful on the one beside it Iterating, instead of factorising

A parameter that is also a price

ξ = coth(Pe) − 1/Pe is the fraction of h/2 that makes a boundary-layer solution exact at every node. On a problem with no layer in it, the error the same scheme commits is ξ times upwinding's — 0.2511 against a ξ of 0.2504, 0.7461 against 0.7448 — so the number that buys the exactness is also the invoice.

 σ14.87σ22.55σ31.08σ40.782σ58.04·10⁻¹⁷σ610⁻¹⁸σ710⁻¹⁸n = 7, and det(A − λB) has degree 4an integer, and a judgementdegree of det(A − λB), exactly4infinite eigenvalues, from the degree3singular values below the cut3largest gap in the spectrum∞a degree cannot be nearly threeand a singular value can be nearly zero Eigenvalues, singular values, rank

The largest gap is inside the null space

The rule recommended for counting a pencil's infinite eigenvalues is to cut at the largest gap in the singular values of B. On integer pencils, with no perturbation anywhere and an exact answer available from the characteristic polynomial, it returns the wrong count on nine of twenty-five — because the singular values that are mathematically zero come back spread over a hundred and forty orders of magnitude, and the largest ratio in the list is between two of them.

 012345678910⁻¹³10⁻¹⁰10⁻⁷10⁻⁴10⁻¹rank budgetresidual the run stalls atb with no structure at allb constant: the answer is a traintwo ladders, one truncationstructured, rank 10.12structured, rank 45.5·10⁻⁶structured, rank 85.9·10⁻¹⁴unstructured, rank 10.97unstructured, rank 80.47the floor is not the truncation'sit is the answer's Iterating, instead of factorising

A run that is over at step five

A conjugate gradient whose every iterate is cut to a rank budget reaches the floor that budget allows at step 5, 36, 42 or 59, and then does nothing for the rest of the run. Four times the iterations move the floor by a factor of 1.8, and past the answer's own rank they move it the wrong way.

 0183654729010812610⁻¹⁶10⁻¹³10⁻¹⁰10⁻⁷10⁻⁴10⁻¹products with Asizeresidual boundtrue errorbounded memorybasis vectors kept8products with A140worst error in the k wanted7.1·10⁻¹⁵the bound is free and the error is notand the basis never grows Eigenvalues, singular values, rank

The same budget, spent five ways

A restarted method has one budget — products with A — and two ways to spend it, in many short cycles or a few long ones. At about a hundred and forty products the answer is the same to a factor of seven whichever split is chosen, and the residual bound the method reports spans ten orders of magnitude across the same five runs.

 00.7853981.57082.356193.1415900.250.50.751frequency θdamping |g(θ)|the oscillatory half →⅓ — the symmetric optimumno convectionwith convectionits imaginary parta modulus, not a valuesmoothing factor at this ω0.71best over every ω0.71the symmetric operator's, at ω = 2/30.33the imaginary part does not depend on ωso no ω removes it Iterating, instead of factorising

A smoother that stops being one

Weighted Jacobi's smoothing factor on the convection–diffusion operator is a function of the cell Péclet number and nothing else — identical to eight digits at five grid sizes at matched Pe. It is 0.3335 at Pe = 0.016, exactly 1/√2 at Pe = 1, and 5.2190 at Pe = 7.8, where the sweep amplifies the modes it exists to remove.

 -213284300.250.50.751rotation of the anisotropy (degrees)factor / couplingusableconvergence factoraxis couplingdiagonal couplingthe standard answer, and the anglefactor at 0°0.19factor at 45°0.8axis ÷ diagonal coupling at 45°2the hierarchy reads the matrixand the matrix lost the direction Iterating, instead of factorising

How much direction there was to lose

At 45° the nine-point stencil hands smoothed aggregation the same wrong hierarchy at every anisotropy — six strong neighbours per interior point, 121 aggregates, the identical partition from ε = 10⁻⁴ to 0.099. The convergence factor that one hierarchy produces runs from 0.802 to 0.581 over the same range.

 10⁻¹100.20.40.60.8strength threshold θresidual reduction per cycleθ = εsemi-coarseningkept whole rowsfull coarseningone parameter, two methodsbest factor above ε0.053best factor below ε0.1the ratio across the switch231×31 anisotropic operatora switch, not a dial Iterating, instead of factorising

The switch does not know which side is better

The strength threshold moves the coarsening from full to semi at θ = ε exactly, at every anisotropy. Which of the two converges faster is a separate question with a separate answer, and it changes sign between ε = 0.33 and ε = 0.34 — where nothing whatever happens to the switch.

 10⁻³10⁻²10⁻¹1024681012size of the negative eigenvalue, −λproducts before the test firesharder to find, and milder8 spectra, n = 50products at the largest λ3products at the smallest10smallest share of λ recovered0.14largest0.34the one that hidesis the one that matters least Iterating, instead of factorising

A proof that does not ask how large the matrix is

Proving a Hessian indefinite costs three matrix–vector products when the negative eigenvalue is 3 and nine to eleven when it is a thousandth, and that pair of numbers barely moves across a fourfold range in n. The factorisation that settles the same question costs a third of n³, which grows by a factor of sixty-four over the same range.

 10³10⁶10⁹10¹²10¹⁵10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³1largest iterate on the way, ‖x‖difference between the two residualsthe bound, linear in ‖x‖what it does, slope 0.46a bound of one, a walk of a halffitted slope0.46the bound's slope1smallest gap measured3.9·10⁻¹⁵largest gap measured1.5·10⁻¹⁰share of the bound, near end9.1·10⁻⁴share of the bound, far end3.4·10⁻⁹the bound is rightand loose by a square root Iterating, instead of factorising

A walk needs a length

The gap between the two residuals grows as the square root of something, and a square root needs a length. Two quantities are candidates — how far the iterates travelled and how many steps were taken — and only a second sweep separates them. Across a fourfold change in size the iteration count goes from 39 to 96 and the gap goes from 5.04·10⁻¹⁵ to 5.33·10⁻¹⁵.

 22.32.62.910⁻⁴10⁻³10⁻²10⁻¹1log₁₀ numbers storeddistance to the dominant eigenvalueArnoldi, linearisedprojected quadraticper number heldstorage, linearised832storage, second-order416Ritz values, linearised26Ritz values, second-order52half the storageand twice the approximations Iterating, instead of factorising

The answer that arrives when the space runs out

A second-order Krylov recurrence holds vectors of length n for a problem with 2n eigenvalues, so it is exact at n steps where the linearised route needs 2n. The machine-precision reading at forty-four vectors on a chain of forty is that exhaustion rather than convergence, and it arrives through a basis whose ‖QᵀQ − I‖ is above one.

 110¹00.30.60.91.2trust-region radius Δshare of the exact model decreasethe exact subproblem7 radii, n = 60share at the smallest radius0.95share at the largest0.3products, at most8radii stopped by the curvature4a few products against an eigendecompositionand most of the decrease Iterating, instead of factorising

The certificate that arrives soonest is worth least

The more negative a Hessian's smallest eigenvalue, the sooner conjugate gradients meets a direction of negative curvature — and the less of the exact trust-region decrease that direction turns out to be worth. At λₘᵢₙ = −10 the step arrives after two products and gets 39.6 per cent; at −10⁻³ the same two products get 89.8, and the whole sweep costs eight.

All themes