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.

Worth reading first: The units the matrix is measured in · The condition number is an amplifier · The exact answer to a nearby problem.

The essay that introduced the componentwise condition number did so as a repair for an embarrassment: κ₂ can be changed by a factor of 10⁶ by measuring one variable in millimetres instead of kilometres, which means it is partly a statement about somebody’s units. The componentwise number is invariant to that rescaling, so it is the honest version of the same quantity.

That reading makes it sound like a refinement — usually a bit smaller, occasionally worth computing. On one matrix in this collection the two numbers are 3.04·10¹³ and 13.25, and only one of them says anything about the error.

The augmented step's two condition numbers, and the one the error obeysκ₂ is a worst case over all perturbations of the same NORM, and a normwise perturbation is allowed to put its whole budget on the smallest entry of a matrix. The componentwise number — Skeel's ‖ |A⁻¹||A||x| ‖ / ‖x‖ — is a worst case over perturbations proportional to the entries, which is what a backward-stable factorisation actually makes. For the augmented interior-point matrix the two are 3.044·10¹³ and 13.25 at μ = 10⁻¹², a ratio of 2.296·10¹², and the measured error is 9.434·10⁻¹⁶ — which is the second number times the unit roundoff and has nothing to do with the first. The componentwise line is still finding its level down to about μ = 10⁻⁴, where the two groups of constraints have not yet separated, and from there it does not move in the fourth digit; where it settles — 13.25 here — is set by how many constraints are active rather than by μ. Eliminating the second block to reach the condensed form destroys the distinction: there the componentwise number is 1.324·10¹⁴, tracking the normwise one, because the large entries have been summed into CᵀDC and are no longer separately identifiable.-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
Fig. 1 One matrix, two condition numbers, thirteen orders apart, and a measured error that obeys the lower line.

The two definitions, side by side

The normwise condition number of a solve is

κ(A) = ‖A‖ · ‖A⁻¹‖

and the bound it supports is: if the matrix is perturbed to A + ΔA with ‖ΔA‖ ≤ ε‖A‖, the solution moves by at most κ(A)·ε relatively. The hypothesis is a bound on the norm of the perturbation and says nothing about where it is.

Skeel’s componentwise number is

cond(A, x) = ‖ |A⁻¹| |A| |x| ‖ / ‖x

and the bound it supports is: if each entry is perturbed by at most ε times its own magnitude — |ΔA_ij| ≤ ε|A_ij| — the solution moves by at most cond(A, x)·ε relatively. The hypothesis is a bound on each entry relative to that entry.

The second is always at most the first times a scaling factor, and it can be very much smaller. When it is, the difference is a fact about the matrix’s scales, and the ratio is a measure of how much of κ₂ is a statement about a perturbation nobody would ever make.

Two perturbation bounds and the error that was measured, on a 8×8 matrix spread over 8 decades of unitsThree curves against the size of an entrywise relative perturbation. The normwise bound κ_∞·ε is a valid bound and sits 4·10⁷ times above the componentwise one cond(A, x)·ε, which is also a bound and is nearly attained by the worst of forty random perturbations at each size.10⁻¹⁴10⁻¹³10⁻¹²10⁻¹¹10⁻¹⁰10⁻⁹10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³1relative size of the entrywise perturbationrelative forward errorκ_∞ · εcond(A,x) · εmeasuredboth bounds holdκ_∞(A)1.9·10⁸cond(A, x)4.8ratio of the bounds4·10⁷both curves above the data are boundsand only one of them is a measurement
Fig. 2 The two bounds in the field’s own earlier essay, on a matrix constructed so that the difference is visible.

The matrix

The constraint field’s fifth essay supplies it. A Newton step of an interior-point method, written in its unreduced form

[ H Cᵀ ] [ C −D⁻¹ ]

where D = Z/S is a diagonal whose entries run to 1/μ on the active constraints and to μ on the inactive ones, so −D⁻¹ has entries spanning 1/μ². As μ falls the diagonal separates by twenty-four orders, which is what makes the matrix ill conditioned in the normwise sense.

The entries are not uniform in scale and they are not uniform by construction: the matrix has a block of ordinary numbers, a block of enormous ones and a block of tiny ones, and which is which is known. That is exactly the situation the two definitions disagree about.

Both forms of an interior-point step: two condition numbers that climb together and two errors that do notOne Newton step of an interior-point method for a quadratic programme with 8 unknowns and 6 constraints, 3 of them active, written twice. The augmented form is [[H, Cᵀ], [C, −D⁻¹]] and the condensed form is H + CᵀDC, and they have the same solution. As the barrier parameter μ falls, the diagonal D separates — z/s runs to 1/μ on the active constraints and to μ on the inactive ones — so both condition numbers climb: 29.04 to 3.044·10¹⁵ for the augmented form and 155.8 to 2.403·10¹⁶ for the condensed one, within a factor of 7.9 of each other. The two forward errors, both measured against the exact rational solution of the system the machine is holding, do not follow: the condensed form's rises to 0.3098 and the augmented form's stays at 1.035·10⁻¹⁵ — fifteen correct digits at a condition number of 3.04·10¹⁵.-14-12-10-8-6-4-2010⁻¹⁷10⁻¹³10⁻⁹10⁻⁵10⁻¹10³10⁷10¹¹10¹⁵log₁₀ μ — the barrier parametercondition number, and relative errorκ₂, condensedκ₂, augmentederror, condensederror, augmentedagainst a BigInt answerκ₂ augmented, μ = 10⁻¹⁴3·10¹⁵its relative error10⁻¹⁵κ₂ condensed2.4·10¹⁶its relative error0.31the same step, written two waysand only one of them is solvable
Fig. 3 The same matrix from the other essay, with the consequence rather than the explanation.

The measurement

Both numbers, at four barrier parameters, on the augmented form and on the condensed form:

μ κ₂ augmented componentwise componentwise, condensed 10⁰ 2.90·10¹ 22.54 1.35·10² 10⁻⁴ 3.02·10⁵ 13.18 1.31·10⁶ 10⁻⁸ 3.04·10⁹ 13.25 1.32·10¹⁰ 10⁻¹² 3.04·10¹³ 13.25 1.32·10¹⁴

For the augmented matrix the componentwise number is 13.25 and does not move while κ₂ climbs thirteen decades. For the condensed matrix the two agree to within a factor of two at every μ.

And the measured error follows the componentwise number: 9.4·10⁻¹⁶ for the augmented form at μ = 10⁻¹², which is 13.25 · u to within four orders, against a κ₂·u of 3.4·10⁻³.

So the two numbers are not a refinement of one another. On this matrix they are different answers to different questions, and one of the questions is the one being asked.

The augmented step's two condition numbers, and the one the error obeysκ₂ is a worst case over all perturbations of the same NORM, and a normwise perturbation is allowed to put its whole budget on the smallest entry of a matrix. The componentwise number — Skeel's ‖ |A⁻¹||A||x| ‖ / ‖x‖ — is a worst case over perturbations proportional to the entries, which is what a backward-stable factorisation actually makes. For the augmented interior-point matrix the two are 2.52·10¹³ and 11.19 at μ = 10⁻¹², a ratio of 2.252·10¹², and the measured error is 6.091·10⁻¹⁶ — which is the second number times the unit roundoff and has nothing to do with the first. The componentwise line is still finding its level down to about μ = 10⁻⁴, where the two groups of constraints have not yet separated, and from there it does not move in the fourth digit; where it settles — 11.19 here — is set by how many constraints are active rather than by μ. Eliminating the second block to reach the condensed form destroys the distinction: there the componentwise number is 2.759·10¹⁴, tracking the normwise one, because the large entries have been summed into CᵀDC and are no longer separately identifiable.-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⁻¹²2.5·10¹³componentwise, same matrix11their ratio2.3·10¹²measured error6.1·10⁻¹⁶every library prints the top lineand the error obeys the third
Fig. 4 At ten constraints, where the flat line sits a little higher — it rises with the size of the problem and not with μ.

What “13.25” is a number about

A flat line at 13.25 across thirteen decades invites the question of where 13.25 comes from, and the answer is that it is about the size of the problem rather than about the barrier.

cond(A, x) = ‖ |A⁻¹||A||x| ‖ / ‖x‖ is at least one for any matrix, and for a matrix with no cancellation in the products it is a small multiple of the number of terms being summed. Here the matrix is 14 × 14, and 13.25 is what a fourteen-term sum of same-signed contributions produces. Raise the constraint count to ten and the matrix is 18 × 18 and the number goes to about seventeen; lower it to four and it falls.

So the line is flat because nothing about the barrier enters it, and its height is a property of the shape. That is worth having explicitly, because it says the good behaviour is not a lucky cancellation that might not survive a different problem: it is what the componentwise number always is for a matrix whose entries do not conspire, and the barrier does not make them conspire.

The one place it is not flat is at μ = 1, where it reads 22.5 rather than 13.2. At that end the two groups of constraints have not separated, so |A⁻¹||A| has contributions of comparable size from more places, and the count of terms that matter is larger. The separation helps the componentwise number, which is the exact reverse of what it does to κ₂.

The augmented step's two condition numbers, and the one the error obeysκ₂ is a worst case over all perturbations of the same NORM, and a normwise perturbation is allowed to put its whole budget on the smallest entry of a matrix. The componentwise number — Skeel's ‖ |A⁻¹||A||x| ‖ / ‖x‖ — is a worst case over perturbations proportional to the entries, which is what a backward-stable factorisation actually makes. For the augmented interior-point matrix the two are 1.063·10¹³ and 11.09 at μ = 10⁻¹², a ratio of 9.585·10¹¹, and the measured error is 8.739·10⁻¹⁶ — which is the second number times the unit roundoff and has nothing to do with the first. The componentwise line is still finding its level down to about μ = 10⁻⁴, where the two groups of constraints have not yet separated, and from there it does not move in the fourth digit; where it settles — 11.09 here — is set by how many constraints are active rather than by μ. Eliminating the second block to reach the condensed form destroys the distinction: there the componentwise number is 1.803·10¹⁴, tracking the normwise one, because the large entries have been summed into CᵀDC and are no longer separately identifiable.-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⁻¹²1.1·10¹³componentwise, same matrix11their ratio9.6·10¹¹measured error8.7·10⁻¹⁶every library prints the top lineand the error obeys the third
Fig. 5 At four constraints, where the flat line sits lower — because it counts terms and there are fewer of them.

Why the perturbation a factorisation makes is componentwise

The reason the smaller number is the relevant one is a fact about backward stability rather than about condition numbers, and it is worth stating precisely because the usual statement is normwise and hides it.

A backward-stable factorisation returns the exact factors of A + ΔA with ‖ΔA‖ ≤ cn·u·‖A‖. That is the normwise statement, and it is what the backward error essays quote. The finer statement, which holds for Gaussian elimination and for symmetric indefinite factorisation with a bounded growth factor, is componentwise: |ΔA_ij| ≤ cn·u·(|L||U|)_ij, and in the absence of growth that is |ΔA_ij| ≲ cn·u·|A_ij|.

Each entry is perturbed by an amount relative to itself. The enormous entries of −D⁻¹ are perturbed by enormous absolute amounts that are tiny relative to themselves; the ordinary entries of H and C are perturbed by ordinary tiny amounts. That is a componentwise perturbation, and the number that bounds its effect is the componentwise one.

A normwise perturbation of the same norm is permitted to put its entire budget on the smallest entry — a change of 10⁻² to an entry of size 10⁻² — which no factorisation would ever make and which is what κ₂ = 3·10¹³ is a worst case over.

Backward error with and without the structure, for Levinson and for elimination, n = 12Both solvers are backward stable on this family by the usual measure: the smallest perturbation of any kind that explains either answer stays below 3.1·10⁻¹⁷ at every ρ, and the two curves are indistinguishable. Insisting the perturbation be a symmetric Toeplitz matrix — the kind of object the problem was posed with — lifts both by orders, and at ρ = 0.999 Levinson's rises by a factor of 2.48·10⁶. Neither number is alarming in absolute terms; what is worth carrying is that the reassuring one is the one that is reported, and the two are not measuring the same thing.10⁻³10⁻²10⁻¹110⁻¹⁹10⁻¹⁷10⁻¹⁵10⁻¹³10⁻¹¹10⁻⁹1 − ρbackward errorLevinson, Toeplitzelimination, ToeplitzLevinson, any kindelimination, any kindat ρ = 0.999Levinson, any perturbation3.1·10⁻¹⁷Levinson, Toeplitz only7.6·10⁻¹¹the ratio between them2.5·10⁶diagonal defect of the first0.97the number that is reportedand the number that is asked about
Fig. 6 The backward errors several solvers actually produce, in the field’s own essay.

And why the elimination destroys it

The condensed form is H + CᵀDC, and its two condition numbers agree. The structure that protected the augmented form is gone, and the reason is arithmetic rather than algebraic.

An entry of CᵀDC is a sum over the constraints: Σₖ C_ki D_k C_kj. Some terms in that sum are of size 1/μ and some are of size μ, and the result is dominated by the first. So an entry of the condensed matrix is a huge number that has a small number hidden inside it, and perturbing the entry by ε times itself perturbs the hidden contribution by ε times the huge number — which is enormous relative to the hidden part.

The information about scale was in the separateness of the entries, and the elimination summed them. What is lost is not accuracy at the moment of summing — the sum is computed to full relative precision — but the ability of any later perturbation to respect the two scales, because after the sum they are no longer distinguishable.

That is a general lesson about eliminations and it is not confined to this family: an elimination is a change of variables, and a change of variables can destroy structure that the original variables had. The tensor field found the same shape when an all-orthogonality measurement was scaled by the wrong thing, and the structured backward error essays found it when a nearby problem turned out not to be a problem of the right kind.

Both forms of an interior-point step: two condition numbers that climb together and two errors that do notOne Newton step of an interior-point method for a quadratic programme with 8 unknowns and 6 constraints, 1 of them active, written twice. The augmented form is [[H, Cᵀ], [C, −D⁻¹]] and the condensed form is H + CᵀDC, and they have the same solution. As the barrier parameter μ falls, the diagonal D separates — z/s runs to 1/μ on the active constraints and to μ on the inactive ones — so both condition numbers climb: 28.23 to 1.439·10¹⁶ for the augmented form and 162.7 to 4.182·10¹⁷ for the condensed one, within a factor of 29.1 of each other. The two forward errors, both measured against the exact rational solution of the system the machine is holding, do not follow: the condensed form's rises to 1 and the augmented form's stays at 1.352·10⁻¹⁵ — fifteen correct digits at a condition number of 1.44·10¹⁶.-14-12-10-8-6-4-2010⁻¹⁷10⁻¹³10⁻⁹10⁻⁵10⁻¹10³10⁷10¹¹10¹⁵log₁₀ μ — the barrier parametercondition number, and relative errorκ₂, condensedκ₂, augmentederror, condensederror, augmentedagainst a BigInt answerκ₂ augmented, μ = 10⁻¹⁴1.4·10¹⁶its relative error1.4·10⁻¹⁵κ₂ condensed4.2·10¹⁷its relative error1the same step, written two waysand only one of them is solvable
Fig. 7 Where that leads, measured: two errors, fourteen orders apart, on a pair of condition numbers a factor of eight apart.

What each number is good for

Neither is the right number in general, and saying which question each answers is more useful than picking a winner.

κ₂ is right for an error in the data. A measured H, a C assembled from noisy geometry, a right-hand side from an instrument: those errors are not proportional to the entries they land on, and against them the matrix really is as sensitive as κ₂ says. A page that concluded “κ₂ is misleading” would be wrong about the case it was invented for.

The componentwise number is right for an error the arithmetic makes, which is the case this site is about. It is also right for an error in data that is known relatively — a measurement quoted to three significant figures perturbs each entry proportionally.

And there is a third case neither covers: a structured perturbation that has to preserve a pattern. The structured backward error of a Toeplitz solve is a different number again, because the perturbation is required to stay Toeplitz. Each of the three is a worst case over a different set, and the sets are nested.

κ(A) and κ of the row-equilibrated matrix, along both families, n = 8Write A = D·X with D diagonal and every row of X of unit norm. Every theorem about relative accuracy is a hypothesis on κ(X), and κ(A) appears in none of them — which is easy to read past and is the whole difference between the two families here. Along the graded family κ(A) climbs from 271 to 7.04·10⁵⁰ and κ(X) is 4.892 at every one of the six matrices — the grading is exactly what the diagonal factor absorbs. Along the uniform family κ(X) climbs to 3.87·10²⁵. That is the number that says which question has an answer, and it is not the number anybody prints.01020304050110¹⁰10²⁰10³⁰10⁴⁰10⁵⁰log₁₀ κ(A)log₁₀ of the condition numberκ(X) = κ(A)graded familyuniform familyone number decides, and it is not κκ(X), graded, at every grading4.9κ(A), graded, at the widest7·10⁵⁰κ(X), uniform, at the widest3.9·10²⁵κ(A) there8·10²⁸κ is a fact about the matrixand the hypothesis is about a factor of it
Fig. 8 The case that made the componentwise number necessary, from this field’s earlier essay: a condition number that a change of units moves and one that it cannot.

The third quantity in the identity

The site’s spine is forward ⪅ condition × backward, and this page has been about replacing one factor. The other factor changes too, and the two changes are the same change seen twice.

The normwise backward error of a computed x is ‖Axb‖ / (‖A‖‖x‖ + ‖b‖), and its componentwise cousin divides each component of the residual by (|A||x| + |b|)ᵢ instead. For the augmented matrix the second is at the rounding level while the first is too, because both are measuring an algorithm that behaved well — the backward error is the easy half here.

What the pairing means is that the identity has to be used consistently: a componentwise condition number multiplies a componentwise backward error, and a normwise one multiplies a normwise one. Mixing them gives a bound that is not a bound in either sense. The reason the normwise product is so far above the measured error on this matrix is not that the bound is weak — it is a perfectly good bound on the wrong quantity — but that the pair being multiplied does not describe what happened.

That is the same failure this field has recorded once before, from the other direction, when a backward error turned out to promise a nearby problem of the wrong kind. The identity is exact; what it is about depends entirely on which set of perturbations both halves are worst cases over.

Distance from the truth and residual, for both methods, with 0% of the noise in the matrixFour bars, medians over 40 seeds. The upper pair is how far each answer is from the coefficients the problem was built from; the lower pair is ‖Ax − b‖ on the problem as given. Ordinary least squares minimises the lower quantity by definition, so its bar is the shorter of the two whatever happens above — and at this share it is the less accurate answer.the upper pair is distance from the truth; the lower pair is ‖Ax − b‖least squares · error0.04439total least squares · error0.1565least squares · ‖Ax − b‖4.734total least squares · ‖Ax − b‖4.971two orderingserror ratio (ls ÷ tls)0.28residual ratio (tls ÷ ls)1.1seeds40no vector makes the residual smallernot even the one the problem was built from
Fig. 9 The pairing done consistently, in the field that owns the structured version.

What it costs to compute

cond(A, x) needs |A⁻¹|, and forming an inverse is the operation this collection has an essay against. So the definition is not directly usable, which is one reason libraries print κ₂ instead.

What is usable is an estimate. The quantity ‖ |A⁻¹||A||x| ‖ is the norm of A⁻¹ applied to a known nonnegative vector, so it is one solve with the already-computed factors per estimate, and the standard condition estimator that this field has already priced applies unchanged. That makes it about as expensive as the normwise estimate every library already computes.

The measurement on this page forms the inverse, because the matrices are fourteen by fourteen and the point is the exact value rather than the practicality. On a real problem the estimate is what would be used, and it inherits the estimator’s own failure mode — a matrix on which it under-reports by an arbitrary factor.

How close Hager's estimate is to the true κ₁, over 200 seeded 8×8 matricesFive bars. The estimate is exactly the true condition number on 81% of the sample and inside ten per cent on 87%; the worst underestimate in the whole sample returns 38% of the truth. The last bar is the matrix built to defeat it, at 7.7%, well below anything the random sample reached.each bar is a percentage — of the sample, or of the true condition numberexactly right80.5%inside 10%87.0%inside a factor of 287.0%worst in the sample, ×10037.7%the constructed matrix, ×1007.7%usually exactexact share0.81worst of the sample0.38the constructed matrix0.077a routine that is right most of the timeand never wrong in the safe direction
Fig. 10 The estimator, and the matrix it can be fooled by.

Which of the two a reader should print

The practical conclusion is not “use the componentwise number”, because on most matrices the two agree and the extra solve is wasted. It is a rule about when they can disagree, and the rule has a shape that is easy to check.

They can disagree only when the entries of the matrix span many orders of magnitude and the spread is structural rather than accidental — when the matrix has blocks or rows or columns of different scale, and those scales came from the problem’s construction rather than from noise. A matrix whose entries are all within a few orders has cond(A, x) within a small factor of κ∞(A) and nothing to learn.

That test is cheap: the ratio of the largest entry to the smallest nonzero one, which any code already knows. If it is under 10⁶ the two numbers agree; if it is 10²⁴, as it is here, they may not, and the componentwise estimate is worth the solve.

The same test says where else on this site the question would have been worth asking. A graded matrix in the relative-accuracy essays has exactly this shape, and its eigenvalues turn out to have relative accuracy nobody expected for exactly this reason; so does a scaled system in the units essay. Three places, one diagnostic.

The worst relative error over the whole spectrum, against how many decades the matrix is graded over, n = 8Each point is one bidiagonal matrix and the worst relative error any of its singular values suffers, measured against a Sturm bisection in exact rationals. The route through BᵀB is as good as anything at 2.1 decades — 3.63·10⁻¹⁴ — and by 11 decades it is at 14.5, which is not an error in the answer, it is the answer. One-sided Jacobi and the zero-shift bidiagonal sweep are flat at about the unit roundoff across the whole range, which is the claim this site's own SVD has been making in a source comment since it was written and had never measured.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
Fig. 11 A graded matrix in the field where the same spread turns out to buy something.

The refusal

The reading to close is the one this essay’s own predecessor invites: the componentwise condition number is a refinement of the normwise one — usually a bit smaller, occasionally worth computing, and never a different answer. Every word of that is a reasonable summary of the earlier essay, and it is what makes the case here worth measuring.

The assertion is fed both numbers for the augmented matrix at μ = 10⁻¹² and required to reject the claim that they agree within six orders. They are 3.04·10¹³ and 13.25, a ratio of 2.3·10¹².

Two further refusals guard the other direction. One is fed the claim that eliminating a block keeps the componentwise conditioning, and required to refuse — the whole second half of the page rests on the elimination being what destroys it. The other is fed the claim that a condition number of 10¹³ costs thirteen digits, and required to refuse: the standard bound is an upper bound and this matrix is twelve orders under it.

The augmented step's two condition numbers, and the one the error obeysκ₂ is a worst case over all perturbations of the same NORM, and a normwise perturbation is allowed to put its whole budget on the smallest entry of a matrix. The componentwise number — Skeel's ‖ |A⁻¹||A||x| ‖ / ‖x‖ — is a worst case over perturbations proportional to the entries, which is what a backward-stable factorisation actually makes. For the augmented interior-point matrix the two are 1.063·10¹³ and 11.09 at μ = 10⁻¹², a ratio of 9.585·10¹¹, and the measured error is 8.739·10⁻¹⁶ — which is the second number times the unit roundoff and has nothing to do with the first. The componentwise line is still finding its level down to about μ = 10⁻⁴, where the two groups of constraints have not yet separated, and from there it does not move in the fourth digit; where it settles — 11.09 here — is set by how many constraints are active rather than by μ. Eliminating the second block to reach the condensed form destroys the distinction: there the componentwise number is 1.803·10¹⁴, tracking the normwise one, because the large entries have been summed into CᵀDC and are no longer separately identifiable.-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⁻¹²1.1·10¹³componentwise, same matrix11their ratio9.6·10¹¹measured error8.7·10⁻¹⁶every library prints the top lineand the error obeys the third
Fig. 12 At four constraints, where the flat line is lower and the separation is the same.

A note on which norm

Both definitions above were written with a norm left unspecified, and the choice is not free.

The componentwise bound is naturally an infinity-norm statement, because it comes from bounding each component of the error separately and then taking the largest. Mixing it with a two-norm condition number — which is what a code that prints κ₂ and then quotes a componentwise bound would be doing — introduces a factor of √n that is harmless at n = 14 and is not at n = 10⁶.

The measurements on this page use the infinity norm throughout for the componentwise number and the two-norm for κ₂, which is the convention each is usually quoted in, and the ratio between them for a matrix of this size is under four. That is well inside the thirteen orders being reported, so nothing here turns on it — but a reader carrying the numbers elsewhere should carry the norms with them.

Two perturbation bounds and the error that was measured, on a 8×8 matrix spread over 8 decades of unitsThree curves against the size of an entrywise relative perturbation. The normwise bound κ_∞·ε is a valid bound and sits 4·10⁷ times above the componentwise one cond(A, x)·ε, which is also a bound and is nearly attained by the worst of forty random perturbations at each size.10⁻¹⁴10⁻¹³10⁻¹²10⁻¹¹10⁻¹⁰10⁻⁹10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³1relative size of the entrywise perturbationrelative forward errorκ_∞ · εcond(A,x) · εmeasuredboth bounds holdκ_∞(A)1.9·10⁸cond(A, x)4.8ratio of the bounds4·10⁷both curves above the data are boundsand only one of them is a measurement
Fig. 13 The componentwise bound at a wider spread of scales, where the choice of norm starts to matter.

At other settings

The augmented step's two condition numbers, and the one the error obeysκ₂ is a worst case over all perturbations of the same NORM, and a normwise perturbation is allowed to put its whole budget on the smallest entry of a matrix. The componentwise number — Skeel's ‖ |A⁻¹||A||x| ‖ / ‖x‖ — is a worst case over perturbations proportional to the entries, which is what a backward-stable factorisation actually makes. For the augmented interior-point matrix the two are 1.822·10¹³ and 15.57 at μ = 10⁻¹², a ratio of 1.17·10¹², and the measured error is 6.124·10⁻¹⁶ — which is the second number times the unit roundoff and has nothing to do with the first. The componentwise line is still finding its level down to about μ = 10⁻⁴, where the two groups of constraints have not yet separated, and from there it does not move in the fourth digit; where it settles — 15.57 here — is set by how many constraints are active rather than by μ. Eliminating the second block to reach the condensed form destroys the distinction: there the componentwise number is 4.19·10¹⁴, tracking the normwise one, because the large entries have been summed into CᵀDC and are no longer separately identifiable.-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⁻¹²1.8·10¹³componentwise, same matrix16their ratio1.2·10¹²measured error6.1·10⁻¹⁶every library prints the top lineand the error obeys the third
Fig. 14 Five, for reading against the hero’s six.
The augmented step's two condition numbers, and the one the error obeysκ₂ is a worst case over all perturbations of the same NORM, and a normwise perturbation is allowed to put its whole budget on the smallest entry of a matrix. The componentwise number — Skeel's ‖ |A⁻¹||A||x| ‖ / ‖x‖ — is a worst case over perturbations proportional to the entries, which is what a backward-stable factorisation actually makes. For the augmented interior-point matrix the two are 2.822·10¹³ and 13.17 at μ = 10⁻¹², a ratio of 2.142·10¹², and the measured error is 1.205·10⁻¹⁶ — which is the second number times the unit roundoff and has nothing to do with the first. The componentwise line is still finding its level down to about μ = 10⁻⁴, where the two groups of constraints have not yet separated, and from there it does not move in the fourth digit; where it settles — 13.17 here — is set by how many constraints are active rather than by μ. Eliminating the second block to reach the condensed form destroys the distinction: there the componentwise number is 1.081·10¹⁴, tracking the normwise one, because the large entries have been summed into CᵀDC and are no longer separately identifiable.-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⁻¹²2.8·10¹³componentwise, same matrix13their ratio2.1·10¹²measured error1.2·10⁻¹⁶every library prints the top lineand the error obeys the third
Fig. 15 And eight, where the componentwise number has risen with the count and κ₂ has not noticed.
Both forms of an interior-point step: two condition numbers that climb together and two errors that do notOne Newton step of an interior-point method for a quadratic programme with 8 unknowns and 6 constraints, 2 of them active, written twice. The augmented form is [[H, Cᵀ], [C, −D⁻¹]] and the condensed form is H + CᵀDC, and they have the same solution. As the barrier parameter μ falls, the diagonal D separates — z/s runs to 1/μ on the active constraints and to μ on the inactive ones — so both condition numbers climb: 28.34 to 1.401·10¹⁶ for the augmented form and 160.6 to 1.673·10¹⁷ for the condensed one, within a factor of 11.9 of each other. The two forward errors, both measured against the exact rational solution of the system the machine is holding, do not follow: the condensed form's rises to 1 and the augmented form's stays at 5.024·10⁻¹⁵ — fifteen correct digits at a condition number of 1.4·10¹⁶.-14-12-10-8-6-4-2010⁻¹⁷10⁻¹³10⁻⁹10⁻⁵10⁻¹10³10⁷10¹¹10¹⁵log₁₀ μ — the barrier parametercondition number, and relative errorκ₂, condensedκ₂, augmentederror, condensederror, augmentedagainst a BigInt answerκ₂ augmented, μ = 10⁻¹⁴1.4·10¹⁶its relative error5·10⁻¹⁵κ₂ condensed1.7·10¹⁷its relative error1the same step, written two waysand only one of them is solvable
Fig. 16 The consequence at two active constraints.
How much a perturbation of the right-hand side is amplified, κ = 10⁶The cumulative distribution of the amplification factor over two hundred random perturbation directions, with the condition number marked as the upper limit.110¹10²10³10⁴10⁵10⁶10⁷00.250.50.751amplification of the input perturbationfraction of directions at or belowκ = 10·10⁵worst found 7.6·10⁵6×6, 200 directionsmedian reaches 0.29 of κ
Fig. 17 The normwise bound doing exactly what it was written to do, on a matrix where it is right.
Backward and forward error against the condition numberA log–log plot over twelve decades of condition number. The backward error is a flat line at ten to the minus sixteen; the forward error rises in proportion to the condition number.110²10⁴10⁶10⁸10¹⁰10¹²10¹⁴10⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²10¹condition number κ(A)relative errorforward errorbackward errorpredicted: κ · u8×8, 20 seeds per κ; dashed is the worstthe problem worsens, not the method
Fig. 18 The identity both numbers are used inside.
The backward error of a hierarchical solve, against the error of the representation it usedSolve A_H x = b exactly and the residual against the matrix that was wanted is b − Ax = (A_H − A)x, so ‖b − Ax‖ ⁄ ‖A‖‖x‖ cannot be anything but the representation's own error. Measured across six accuracies spanning ten decades, the two track at a slope of 0.994 and sit a constant 6.7× apart, which is the difference between a norm of a matrix and a norm of that matrix applied to one vector. The reading is the one this field is built on: the accuracy knob is not an accuracy, it is a backward error chosen in advance. Everywhere else on this site a backward error is something an algorithm produced and somebody then measured; here it is a line in the program, and its size is known before the solve starts.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
Fig. 19 The backward error, which is the quantity being multiplied.
Four pivot rules on [[ε, 1], [1, 1]] with ε = 10⁻¹⁶, before and after its first row is scaled by 1/εEight bars of forward error. On the matrix as given, no pivoting loses the answer entirely and all three pivot rules return it to rounding. Multiply the first row by 1/ε — which does not change the solution — and partial pivoting makes no interchange at all and returns the same wrong answer that no pivoting does. Scaled partial pivoting and complete pivoting are unmoved.forward error, relative to a solution of exactly (1, 1)no pivoting · as given1.2 0 interchangesno pivoting · rows scaled1 0 interchangespartial · as given1.1·10⁻¹⁶ 1 interchangepartial · rows scaled1 0 interchangesscaled partial · as given1.1·10⁻¹⁶ 1 interchangescaled partial · rows scaled2.2·10⁻¹⁶ 1 interchangecomplete · as given2.2·10⁻¹⁶ 1 interchangecomplete · rows scaled2.2·10⁻¹⁶ 1 interchangethe same problem twicepartial, as given1.1·10⁻¹⁶partial, rows scaled1its relative residual5·10⁻¹⁷complete, rows scaled2.2·10⁻¹⁶the two systems have the same solutionand one pivot rule cannot see it
Fig. 20 What a change of units does to a pivot choice, from the elimination field.
One 40×40 matrix described four ways, with its condition number and its backward error at eachThe matrix is ρ^|i−j| at ρ = 0.9999, where κ = 7.977·10⁵. Described as n² entries its condition number is that; as 79 constant diagonals it is 3.343·10⁵; as the 40 a symmetric Toeplitz matrix has, 3.218·10⁵; and as the one number ρ that the matrix actually holds, 613.2. The backward error of the same computed solution runs the other way — 1.56·10⁻¹⁷, 7.49·10⁻¹⁵, 1.5·10⁻¹¹ — and at the last rung there is none: no ρ whatever has the computed answer as its exact solution, and the closest one leaves 100 per cent of the residual unexplained. The drop that matters is the last one, and it is not a linear-algebra structure at all.110¹10²10³10⁻¹⁸10⁻¹⁴10⁻¹⁰10⁻⁶10⁻²10²10⁶numbers that describe the matrixcondition number, and backward errordensetoeplitzsymmetricρ aloneno such problemcondition numberbackward errorone matrix, four descriptionsκ, all n² entries8·10⁵as one number, ρ613backward error, unconstrained1.6·10⁻¹⁷as a symmetric Toeplitz matrix1.5·10⁻¹¹fewer numbers, better conditionedand no nearby problem left
Fig. 21 The ladder of increasingly restricted perturbations, from the structure field.
The same three schemes on a problem with no layer in itWorst nodal error against the grid size, both axes logarithmic, for a manufactured smooth solution on the identical operator at ε = 0.005. Adding no diffusion gives 0.0016, 4·10⁻⁴, 10·10⁻⁵, falling by four at each refinement. The tuned diffusion gives 0.067, 0.022, 0.0061 — 42 times worse at the coarsest grid, and falling more slowly.10²10⁻⁴10⁻³10⁻²10⁻¹grid points nworst nodal errorupwindtunedcentralthe same tuning, another problemtuned ÷ central at n = 3142tuned ÷ central at n = 12761central's error at the finest grid10·10⁻⁵exact on the problem it was derived fromand harmful on the one beside it
Fig. 22 And the case where the nearby problem a backward error promises is not a problem of the right kind at all.
Backward and forward error of two routes to x, on 30×30 systems across eight decades of κBoth routes start from the same LU factorisation. The LU solve's backward error is 2.2·10⁻¹⁷ at every conditioning — a flat line at the unit roundoff. Multiplying by the explicitly formed inverse gives 4.5·10⁻⁵ at κ = 10¹⁴, a slope of 0.94 against κ. The two forward errors, drawn above them, are 2.8·10⁻⁴ and 0.015 — within a factor of 54, which is why the difference between the two methods is invisible to anyone measuring how wrong the answer is.10²10⁴10⁶10⁸10¹⁰10¹²10¹⁴10⁻¹⁸10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³1κ₂(A)relative errorforward, A⁻¹bforward, A\bbackward, A⁻¹bbackward, A\bat κ = 10¹⁴η, LU solve2.2·10⁻¹⁷η, via the inverse4.5·10⁻⁵forward, LU solve2.8·10⁻⁴forward, via the inverse0.015one factorisation, two ways to use itand one of them forfeits the backward error
Fig. 23 Why the definition’s |A⁻¹| is not formed in practice.
500 random symmetric orderings of a regularised saddle-point matrix, and what each factorisation reproducesA quasi-definite matrix — [[H + δI, Aᵀ], [A, −δI]] with δ = 10⁻⁶ — has an LDLᵀ factorisation with diagonal D for every symmetric permutation, so the ordering can be chosen for fill with no numerical veto at all. That is the count: 500 of 500 orderings factorise here, against 69.0 per cent for the same matrix with the zero block left where the problem put it. Each bar counts the orderings whose factorisation left a relative residual ‖PKPᵀ − LDLᵀ‖/‖K‖ in that decade, with the tallest holding 136 of the 500. They span 1.56·10⁻¹⁶ to 1.02·10⁻⁷ — 9 orders — and the growth factor across them runs from 1 to 6.43·10⁵ — which is 0.643/δ, an inverse law that holds at every δ this figure is drawn at. Existence is not stability, and the theorem says only the first.-18-15-12-9-6-300log₁₀ ‖PKPᵀ − LDLᵀ‖ / ‖K‖share of orderingsbestworstexistence and stabilityfactorise1unregularised0.69worst growth6.4·10⁵growth × δ0.64the ordering is free to chooseand not free of consequence
Fig. 24 And a factorisation whose growth factor is what the componentwise bound depends on.
The augmented step's two condition numbers, and the one the error obeysκ₂ is a worst case over all perturbations of the same NORM, and a normwise perturbation is allowed to put its whole budget on the smallest entry of a matrix. The componentwise number — Skeel's ‖ |A⁻¹||A||x| ‖ / ‖x‖ — is a worst case over perturbations proportional to the entries, which is what a backward-stable factorisation actually makes. For the augmented interior-point matrix the two are 1.401·10¹⁴ and 63.05 at μ = 10⁻¹², a ratio of 2.222·10¹², and the measured error is 3.676·10⁻¹⁵ — which is the second number times the unit roundoff and has nothing to do with the first. The componentwise line is still finding its level down to about μ = 10⁻⁴, where the two groups of constraints have not yet separated, and from there it does not move in the fourth digit; where it settles — 63.05 here — is set by how many constraints are active rather than by μ. Eliminating the second block to reach the condensed form destroys the distinction: there the componentwise number is 5.214·10¹⁴, tracking the normwise one, because the large entries have been summed into CᵀDC and are no longer separately identifiable.-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⁻¹²1.4·10¹⁴componentwise, same matrix63their ratio2.2·10¹²measured error3.7·10⁻¹⁵every library prints the top lineand the error obeys the third
Fig. 25 With two active constraints, where the componentwise line settles at 63.05 rather than 13.25 — the level is set by the active set and the flatness is not.
Both forms of an interior-point step: two condition numbers that climb together and two errors that do notOne Newton step of an interior-point method for a quadratic programme with 8 unknowns and 6 constraints, 4 of them active, written twice. The augmented form is [[H, Cᵀ], [C, −D⁻¹]] and the condensed form is H + CᵀDC, and they have the same solution. As the barrier parameter μ falls, the diagonal D separates — z/s runs to 1/μ on the active constraints and to μ on the inactive ones — so both condition numbers climb: 31.22 to 2.193·10¹⁵ for the augmented form and 102.3 to 2.794·10¹⁶ for the condensed one, within a factor of 12.7 of each other. The two forward errors, both measured against the exact rational solution of the system the machine is holding, do not follow: the condensed form's rises to 0.06652 and the augmented form's stays at 3.49·10⁻¹⁵ — fifteen correct digits at a condition number of 2.19·10¹⁵.-14-12-10-8-6-4-2010⁻¹⁷10⁻¹³10⁻⁹10⁻⁵10⁻¹10³10⁷10¹¹10¹⁵log₁₀ μ — the barrier parametercondition number, and relative errorκ₂, condensedκ₂, augmentederror, condensederror, augmentedagainst a BigInt answerκ₂ augmented, μ = 10⁻¹⁴2.2·10¹⁵its relative error3.5·10⁻¹⁵κ₂ condensed2.8·10¹⁶its relative error0.067the same step, written two waysand only one of them is solvable
Fig. 26 Four active, where the ill-conditioning lives in four directions.
What a factorisation survives across one barrier step, against the step's reduction factorFrom one interior-point iteration to the next, H, C and the sparsity pattern are identical and exactly 5 entries of the matrix change — the diagonal of the (2, 2) block, and nothing off it. That is the sparsity pattern reused for free forever. The factorisation is a different question, and the answer is a threshold: carrying the decomposition of K(μ) to K(σμ) and cleaning up with iterative refinement holds at the rounding level while σ is within a per cent of one, and by σ = 0.9 it buys one step. At the schedules an interior-point method actually uses — σ between 0.5 and 0.1 — the first reuse is already at 0.522 and 1.48·10¹¹. The reason is in the entry count: the 5 entries that moved are the 5 that dominate the matrix, and they moved by a factor of 1/σ, so the relative change in K across one step at σ = 0.1 is 9 rather than the 10⁻³ the collection's reuse essays are about.-1-0.75-0.5-0.25010⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²10¹10⁴log₁₀ σ — the barrier's reduction factorresidual after one reused stepconvergedthe pattern free, the factors notentries moved5off-diagonal0survived at σ = 0.996survived at σ = 0.106 steps0 stepsthe few entries that movedare the ones that dominate
Fig. 27 The same matrix’s other awkward property, from the sequence field.
What each way of eliminating a constraint inherits, over five decades of κ(A)The same system solved twice, at 8 unknowns and 3 constraints with κ(H) = 10⁶. The range-space method forms S = AH⁻¹Aᵀ and inherits κ(S), which rises from 1700 to 1.118·10¹² — the square of κ(A), for the reason the normal equations square it. The null-space method solves with the reduced Hessian ZᵀHZ, whose condition number is 7646 at the start of the sweep and 7646 at the end: it does not contain κ(A) at all. The two forward errors, measured against a solution computed in BigInt rationals, follow their own condition numbers: 2.467·10⁻⁶ against 4.419·10⁻¹² at the far end. Both methods are correct and one of them is usable.01234510⁻¹⁷10⁻¹³10⁻⁹10⁻⁵10⁻¹10³10⁷10¹¹log₁₀ κ(A)condition number, and relative errorκ(S)κ(ZᵀHZ)range-space errornull-space erroragainst a BigInt answerκ(S) at κ(A) = 10⁵1.1·10¹²κ(ZᵀHZ), all stops7646range-space forward error2.5·10⁻⁶null-space forward error4.4·10⁻¹²both are the same algebraand only one squares
Fig. 28 And an elimination whose damage the normwise number does report correctly.
The pivots Cholesky computes on a saddle-point matrix, and the one it stops at — row 11 of 14Every pivot Cholesky forms on K = [[H, Aᵀ], [A, 0]], in order, with the 10 taken inside H drawn as positive and the 11-th, which is not, drawn as the failure. It is not a near miss and it is not a rounding accident: the leading n × n block is positive definite, so the first 10 pivots are positive at any conditioning, and the next one is 0 − a₁ᵀH⁻¹a₁ where a₁ᵀ is the first row of A. That equals −(AH⁻¹Aᵀ)₁₁ = -328.56, and the pivot the routine actually reports is -328.56 — the same number to 13 digits. The failure has an address and a value, and both were available before the factorisation started.02468101210⁻⁵10⁻³10⁻¹10¹10³pivot, in order|pivot|negative, at row 11|S₁₁| from AH⁻¹Aᵀa failure with an addressrows before it10the pivot-329−(AH⁻¹Aᵀ)₁₁-329disagreement8.8·10⁻¹⁴it does not fail somewhereit fails at the constraint
Fig. 29 The matrix’s unregularised relative, from the constraint field.
How nearly the constraint is satisfied, when a sketch is allowed to see it and when it is notSketching reduces a least-squares problem to 8 rows and keeps its minimiser within (1 + ε) of the original's. That is a statement about a norm, and an equality constraint is a statement that a quantity is zero, which no multiplicative distortion preserves. Kept out of the sketch and imposed exactly, the constraint holds to 1.68·10⁻¹⁶ at every weight. Written as a weight and sketched along with everything else, it holds at 1.747·10⁻⁸ at τ = 10⁸ — and the slope is the finding: the violation falls as 1/τ, one power, where the unsketched weighted problem's falls as 1/τ². The sketch costs exactly half the decades the weight was buying, and buying them back means squaring τ, which the previous figure says the arithmetic stops allowing.2468101210⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²log₁₀ τ‖Bx − d‖ / ‖d‖sketched with the objectiveweighted, not sketchedkept out of the sketchone power instead of twokept out — feasibility1.7·10⁻¹⁶sketched at τ = 10⁸1.7·10⁻⁸unsketched at τ = 10⁸1.7·10⁻¹⁶objective ÷ optimum1.2a sketch preserves a normand a constraint is not one
Fig. 30 And a third worst-case-over-a-set, from the randomised field.

What links here

Computed from the collection, not written here: the essays that point at this one.

Shares its objects with

Essays that name at least two of the same things, and that neither author linked.

Named objects

A flat tag is an object no other essay names yet.

Backward errorComponentwise condition numberCondition numberExact ground truthForward errorInterior point methodSaddle-point systemsScalingSkeel condition numberStructured backward error