Numerical Linear Algebra
where the algebra stops being the arithmetic
Concept

Lyapunov equation — where it appears

A matrix equation in which the unknown appears with the coefficient on one side and its transpose on the other. It is the Sylvester equation with one matrix repeated, so a Bartels–Stewart solver answers it with no new algorithm.

Named by 13 essays across one field — each of them below, with the objects they name alongside it.

123456710⁻¹¹10⁻⁹10⁻⁷10⁻⁵10⁻³10⁻¹order r kepterror and bounds2Σσ, and the errorσᵣ₊₁a bound that is an equalityorders7bound ÷ error, worst1spread over the sweep1error ÷ σᵣ₊₁, worst2.1computed before the modeland attained by it

The bound that is known in advance

Almost every error in numerical linear algebra is measured after the fact. Balanced truncation has one that is computable before the reduced model exists, in a norm of a function rather than of a residual — and on ordinary problems it is not an upper bound that is loose. It is attained.

reduction · Balanced truncation 1357911131510⁻¹⁸10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³1index kσₖ ÷ σ₁σ₁√utwo factorstheir productthe normal equations, againroutes agree to k =6product floors at9.7·10⁻¹⁰σ₁√u2.3·10⁻⁹κ(P)κ(Q)1.8·10³⁶√ of it1.3·10¹⁸do not form the productthe σ below the line are the bound

The product nobody had to form

The Hankel singular values are the square roots of the eigenvalues of PQ. Form that product and half of them stop existing, at a floor one number predicts — and the fix is the one the least-squares field has had since its first essay, arriving in a place with no least-squares problem in it.

reduction · Balanced truncation 01234567891011121310⁻¹110¹10²10³iterationinterpolation point σa fixed point that is a conditionorder4iterations13σ against −λ(Aᵣ)2.8·10⁻¹²H₂ error4.3·10⁻⁵it stopped movingand the condition holds there

Interpolating at the model’s own poles

One choice of interpolation points is not arbitrary — the mirrored poles of the model about to be built. It is a fixed point rather than a guess, and when it is reached it beats a method costing O(n³) — by 0.4 per cent, which is the honest size of the whole contest.

reduction · Moment matching 1234567891010⁻¹⁹10⁻¹⁶10⁻¹³10⁻¹⁰10⁻⁷10⁻⁴10⁻¹kλₖ₊₁ ÷ λ₁, and the boundZₖ²the Gramianpredicted from two numbersstates30κ of the spectrum389λ11 ÷ λ₁2.5·10⁻⁸the bound there5.2·10⁻⁴the cliff everything rests onand the reason for it

Why a Gramian can be truncated at all

Every method in this field rests on one fact nobody states the reason for — the eigenvalues of a Gramian fall off a cliff. The equation defining it has a rank-one right-hand side and no low-rank structure anywhere — and the answer's decay is a rational approximation problem with a closed-form rate.

reduction · Gramian decay 10¹10²10³10⁻⁹10⁻⁷10⁻⁵10⁻³10⁻¹λ, over the spectrum of −A|r(λ)|geometric maxequally spacedgeometricone line of codeshifts6κ of the spectrum389geometric max0.12equally spaced max0.89the factor between7.5the same k solvesand one choice of where

Where to put the poles of a rational function

Three times in one field the same question arrives from different directions — ADI shifts, rational approximation of a square root, the decay of a Gramian — and it has one answer. Cluster them geometrically towards wherever the function is difficult, and the alternative that looks reasonable costs orders.

reduction · Gramian decay condition numbers, on a logarithmic scaleκ(P)5.66·10¹⁸κ(Q)5.71·10¹⁸√(κ(P)κ(Q)) — the SVD route5.68·10¹⁸κ(P)κ(Q) — the product route3.23·10³⁷1/u4.5·10¹⁵the mean, or the productand only one of them fits

A condition number that is not the model's

κ(P)κ(Q) is quoted as the reason one route to the Hankel singular values fails, and it is. It is also read as a measure of how reducible a model is, and it is not — the norm of the Gramian does not move at all as the McMillan degree runs from four to fourteen, and the condition number wanders over a factor of thirty-five with no trend.

reduction · Gramian conditioning 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

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.

reduction · Gramian decay ‖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

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.

reduction · Balanced truncation edge of the numerical rangespectrum, -3.49zero37.9-0.986-3.496.59unstable one-sided models, of nineteen placements051015-3-2-1.5-1-0.500.511.522.53cond 1.4·10⁶cond 112t, the rescaling exponent — the k-th state is multiplied by γ to the power tk; t = 1 symmetrisesorder 1order 2order 3order 6the transfer function is the same at every tonly the inner product changes

A certificate written in coordinates

Every one-sided reduced model of the convection–diffusion system was stable, because its poles cannot leave the numerical range and the range sat in the left half plane. Rescale the state by a diagonal matrix and nothing about the system changes — not its eigenvalues, not its transfer function, not one pole of any two-sided model — but the numerical range moves. At Péclet 40 a rescaling with condition number 260 pushes its edge to +6.6 and ten placements of nineteen return unstable one-sided models. Projecting in the inner product of a Lyapunov solution puts the certificate back in any coordinates.

reduction · Reduced stability Péclet 100, orders 1, 2, 3, 6nodal edge at 1,000, outflow8.9weighted edge at 1,000-0.097unstable nodal models7110¹10²10³-20246810largest cell over smallestright edge of the numerical rangenodal, small cells at the outflownodal, small cells at the inflowcell-size weightedabove the dashed line the certificate is goneit is a sufficient condition, not a failure

The inner product the mesh already computed

A Galerkin reduced model is certified stable when the operator's numerical range sits in the left half plane, and the certificate belongs to the coordinates. The coordinates that broke it before were a constructed rescaling. Grade a convection–diffusion mesh towards its outflow boundary layer — the grading anyone resolving the layer would choose — and nodal values break it too: the range's edge is past zero at a ratio of ten and reaches +8.9 at a thousand, and at Péclet 100 seven reduced models come back unstable. Weight the projection by the cell sizes the discretisation already computed and every model at every grading is stable, with the range's edge back at −0.098.

reduction · Reduced stability ζ = 0.03largest over-statement4.8removing one state1123456789101112345order keptbound ÷ measured errordashed: the bound attaineda mode that rings is not a heat mode

A mode that rings is counted twice

Balanced truncation's bound, twice the sum of the discarded Hankel singular values, was attained at every order on the heat model — a formula for the error rather than a bound on it. On six lightly damped oscillators it is loose at every order but the last, by factors between 1.3 and 5, at every damping from 0.7 to 0.01. The reason is visible in the singular values: a mode that rings contributes two of them, nearly equal, and removing the mode costs about what removing one would. The last mode's bound over its error is 2.02, 2.00 as the damping falls to 0.03 and 0.01. The equality was a property of modes with one singular value each.

reduction · Balanced truncation the pair's two σ added: distance from oneworst, damping to 0.050.053worst at a resonance, damping 0.10.160.40.60.811.21.41.61.8damping ratioreading ÷ measured error0.0050.010.030.10.3order 2order 4order 6order 8order 10shaded: within ten per centa pair costs its two σ, until the modes overlap

Where a ringing mode stops being a pair

A lightly damped mode puts two nearly equal Hankel singular values into balanced truncation, and the reduction error was read as twice the larger of the first pair cut — right to two per cent at damping 0.01, wrong by half at 0.1. Read instead as the pair's two values added, the error is right to within six per cent at every order up to damping 0.05, and the pair's split turns out to be twice the damping. What it cannot do is say when the reading fails. The split predicts nothing: the widest split at damping 0.1 has the best reading. The failures come from somewhere the pair does not look — the other dropped modes' response under the peak, and a cut between two singular values that agree to under one per cent, which left a reduced model with a pole a thousand times closer to the axis than any pole of the system it came from.

reduction · Balanced truncation the censusreductions1452ties431doubled, all ties1010⁻³10⁻²10⁻¹110¹10²0.511.52gap at the cut ÷ smaller neighbouring gaperror ÷ twice the first dropped σdoubled abovetiesshaded: gap under a fifth of its neighboursevery doubled error is at a tie

A tie is a gap smaller than its neighbours

A balanced truncation that cut between two Hankel singular values agreeing to under one per cent left a reduced model with a pole near the axis and twice the error the σ promised, and the warning drawn from it was a rule: never cut inside a tie. Over 1,452 reductions of twelve oscillators, all ten doubled errors are at ties — but three of them have gaps of about two per cent, which a one-per-cent rule lets through, and what they share is a gap under a fifth of the gaps beside it. Moving every such cut one order down never doubles and halves each doubled error, for a median cost of one per cent. It is not free everywhere: at damping 0.5 it costs up to two and a half times. And the near-integrator was never the harm: 32 reductions have one and 22 of them read their σ correctly. The cheapest test of all is one evaluation of the error at zero frequency, which picks out the ten doubled cuts and nothing else.

reduction · Balanced truncation

Named alongside it

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

GramianBalanced truncationHankel singular valuesTransfer functionA-priori boundCondition numberSingular valuesReduced stabilityADI iterationHinf normLow-rank approximationMcMillan degree

All concepts