Feasible and wrong

Constraints arrive from different places. A conservation law, a calibration, a normalisation, a boundary condition — and a program that assembles them does not usually check whether any two of them say the same thing.

Two that say the same thing are not an error — unlike a repeated row of a fit, which influence is decided before the data shows changes the answer by changing how the influence is shared. The problem is still well posed, the feasible set is the same set, and the answer is the same answer; a repeated equation constrains nothing new. What changes is the arithmetic, and it changes by as much as the repetition is close.

As ε falls from 1 to 10⁻¹², κ(B) rises from 8.19 to 4.64·10¹². The null-space route’s error rises with it, from 2.96·10⁻¹⁵ to 1.16·10⁻⁴. The saddle-point route reaches 8.22·10⁻³ and the method of weighting 2.61 — no correct digit in either.

And the flat line along the bottom is the feasibility: how nearly each answer satisfies the constraints it was given. It is 10⁻¹⁵ at every point, including the ones with no digits.

What a caller would check, and why it says nothing #

The natural check on a constrained solve is that the constraints hold. It is cheap — one matrix-vector product — it needs no reference answer, and it is the property the whole construction exists to enforce.

It is also satisfied by construction, whatever else has gone wrong.

The null-space route builds its answer as $Q_1 y_1 + Q_2 y_2$ where $y_1$ comes from a triangular solve with the constraint’s own R factor and $y_2$ lives in the null space. The second half changes nothing about feasibility — $Q_2$ spans the null space, so $B(Q_2 y_2) = 0$ exactly, up to the orthogonality of a Householder Q. So the feasibility measures the first half only, and the first half is a small triangular solve that a well-implemented routine gets right to the rounding level whether or not the answer means anything.

The saddle-point route enforces the constraint as m equations of the system it eliminates, and an elimination with pivoting satisfies the equations it was given to the working precision even when the system is nearly singular. That is what backward stability is: the answer solves a nearby system exactly, and nearby leaves the residual small.

So feasibility is a backward-error quantity and the thing that has gone wrong is a forward error. This collection has the standing statement of that distinction in a small residual is not a small error, and the constrained case is the sharper instance: there the residual is small and the error is large because the matrix is ill-conditioned, and here the residual is small because it cannot be anything else.

What each route actually loses #

The three curves separate, and the separation is the same ordering the previous rungs on this ladder found for a different reason.

The null-space route loses κ(B) times the rounding. At the right-hand end that is κ(B)·u = 1.03·10⁻⁴ and the measured error is 1.16·10⁻⁴ — the same number to within thirteen per cent, and the agreement holds along the whole sweep: at κ(B) = 4.64·10⁴, 4.64·10⁶, 4.64·10⁸ and 4.64·10¹⁰ the errors are 8.59·10⁻¹³, 5.87·10⁻¹¹, 2.47·10⁻⁸ and 1.04·10⁻⁶ against κ(B)·u of 1.03·10⁻¹², 1.03·10⁻¹⁰, 1.03·10⁻⁸ and 1.03·10⁻⁶. Four decades of agreement, which is what makes it a law rather than a coincidence at one point. Its only exposure to B is the triangular solve with R, and a triangular solve loses the conditioning of its matrix once.

The saddle-point route loses more. It reaches 1.21·10⁻³ at κ(B) = 4.6·10⁸, where the null-space route is still at 2.47·10⁻⁸ — five orders apart — because the constraint block enters its matrix as a border and the elimination sees the whole (n + m) × (n + m) object, whose conditioning is not κ(B) but something worse.

The three curves also say what the null-space route is not paying, which is the more useful half. κ(B)·u and not κ(B)²·u: at $\kappa(B) = 4.64 \cdot 10^8$ the square would be $2.2$ and the measured error is 2.47·10⁻⁸, so a reader who expected a constrained solve to square its constraint block the way the normal equations square a fit would be wrong by eight orders. The squaring happens in the route that forms a cross-product and nowhere else.

And the weighting method loses most. It is at 4.81·10⁻² where the null-space route is at 2.47·10⁻⁸. Its weighted rows are τB, so a nearly dependent B gives a weighted block with two nearly identical rows of enormous norm, and the row-ordering effect that page measures is amplified by exactly the factor the dependence introduces.

So the ordering is the same as the previous rung’s — null space, then saddle point, then weighting — and the quantity driving it is different. There it was κ(A) and the forming of a cross-product; here A is well conditioned throughout and it is κ(B) and where the constraint block enters.

The two conditionings do not mix #

The figure’s slider raises κ(A) from 10 to 10⁶ and the picture barely changes shape.

That is worth stating because a reader’s instinct is that the two difficulties should compound. They do not, on these problems: the error is roughly the larger of κ(A)·u and κ(B)·u rather than their product, because the two enter different parts of the computation. κ(B) enters the triangular solve that places the answer in the feasible set; κ(A) enters the least-squares problem solved inside it. Neither multiplies the other.

Which has a practical consequence worth more than the measurement. A caller who sees a bad answer and checks κ(A) will find nothing wrong with it — it is 10² on the hero’s sweep, and every curve on that figure is the constraint block’s doing. The diagnosis requires looking at a matrix the caller may not think of as a matrix at all: a list of side conditions, assembled from different parts of a model, which happens to have a rank.

The pair is the clearest statement of the last section. Raising κ(A) by five orders lifts the curves where κ(B) is small and does nothing where κ(B) is large, because the error is the larger of the two contributions rather than their product — so the left of each picture is the fit’s difficulty and the right is the constraints’, and the crossover between them is visible as a bend.

Where a repeated constraint comes from #

Worth a paragraph, because the construction sounds artificial and the situation is not.

A model assembled from parts has constraints from each part. A conservation law at a node and a normalisation over the same set of variables can be the same equation in different units. A boundary condition applied twice to a shared edge is an exact repetition. A calibration constraint and a physical one can agree to five digits and differ in the sixth, which is ε = 10⁻⁵ in the sweep above and an error of 10⁻¹¹ in the answer — invisible, and present.

And the near case is more likely than the exact one, because an exact repetition is usually caught: a rank-revealing factorisation of B finds it, the R factor has a zero on its diagonal, and a routine that checks for that refuses. A near repetition passes every check and costs the digits, which is the rank decision arriving in a place nobody thinks to make one. The threshold a rank check would need is not available either — there is no gap in B’s singular values to cut at, only a continuum, and the measurement above is that continuum.

What to compute instead #

Two numbers, and both come out of the factorisation the null-space route already performs.

κ(B), from the R factor of $B^{\mathsf T}$. The ratio of the largest to the smallest diagonal entry of R is an estimate of it, it costs nothing, and the measurement above says the error is that number times the rounding. A caller told the constraints are conditioned at 10¹² knows immediately that four digits of the answer are gone, and knows it before the fit is solved.

And κ(AQ₂), which the previous rung is about. Between the two of them every failure mode this ladder has measured is covered: a constraint that repeats another shows in the first, a fit whose difficulty the constraint did not remove shows in the second, and an ordinary problem shows in neither.

What is not covered by either is the assembly loss of the saddle-point route, which is a property of the formulation rather than of the problem — so the three numbers a constrained solve should report are κ(B), κ(AQ₂), and which route was used.

What the answer is being compared against #

Every error on this page is a distance from an answer, and a page whose finding is the routines are wrong has to say where the right answer came from.

It is the optimality conditions of the problem as data: the cross-product AᵀA and the vector Aᵀb accumulated in exact rationals from the stored doubles of A and b, the saddle-point system assembled from those exact values, and the whole thing solved in BigInt. Nothing is rounded between the input and the answer, which is the reference the first of these three measurements had to introduce in order to see the assembly loss at all.

That reference is available at every ε in the sweep, including the smallest, because a nearly dependent constraint block is not a singular one: its exact rank is three, the exact solve succeeds, and the answer it returns is a genuine and unique solution of the constrained problem. The routines’ answers are wrong about a well-posed question, and the exactness of the reference is what makes wrong a measurement rather than a comparison between two approximations.

Why the ordering of the routes is the same twice #

The two sweeps agree about which route to use and disagree about why, and the agreement is worth one paragraph because it is the only general recommendation this ladder produces.

When the difficulty is in A, the null-space route wins because it never forms a cross-product and because the problem it solves has already had the constraint removed from it. When the difficulty is in B, it wins because its only exposure to the constraint block is one triangular solve, while the saddle-point route carries the whole block inside a matrix it eliminates and the weighting method multiplies the block by τ before doing anything else.

Two different mechanisms, the same ordering. That is a stronger recommendation than either sweep on its own, and it is the kind a measurement can give and a derivation cannot: the derivations of all three routes are exact, they differ in no algebraic respect, and what separates them is where each one lets a conditioning enter.

The one thing that would overturn it is a problem where the constraint block is large — many constraints, few unknowns left — since the null-space route’s advantage comes from working in a small space and factoring a small $B^{\mathsf T}$. That regime is not drawn here and is the obvious limit of the recommendation.

The refusals #

The claim worth breaking is the one a careful caller makes: the answer satisfies its constraints to the rounding level, so the routine did its job. The assertion that a feasibility below 10⁻¹² implies a correct answer is fed the pair at ε = 10⁻¹² — feasibility 1.03·10⁻¹⁵, error 1.16·10⁻⁴ — and required to fail.

The second refusal is a control on the sweep itself. Fed ε = 1, where the third constraint is a perfectly ordinary independent one, the claim that a third constraint always costs digits has to fail: the null-space route returns 2.96·10⁻¹⁵ there, which is the rounding level, and it is the same routine on the same problem shape.

And the third refuses the reading that the problem has become ill posed. Fed every ε in the sweep, the claim that the constrained problem has no solution has to fail — it has one at every point, the exact reference computes it, and the routines’ answers are wrong about a question that has a perfectly good answer.

Three ways a constrained fit goes wrong #

This ladder has now measured three, and they want three different diagnostics. Setting them beside each other is the practical summary of the whole anchor.

The formulation squares the problem. The optimality conditions contain AᵀA, and assembling them loses the digits before any solve. The diagnostic is knowing which route was used; the repair is a route that does not form a cross-product.

The constraint did not remove the difficulty. κ(A) is the same in both cases and the answers differ by twelve orders. The diagnostic is κ(AQ₂); the repair is nothing — it is a property of the problem, and what the number buys is knowing how much to trust the answer.

The constraints nearly repeat each other, which is this page. The diagnostic is κ(B); the repair is to notice that two side conditions are saying the same thing, which is a modelling question rather than a numerical one.

All three are invisible in the residual and in the feasibility, which are the two quantities a constrained routine returns. That is the shape this collection keeps finding, stated once more: the outputs of a computation are what was asked for, and what would say whether the answer means anything is a quantity nobody asked for and every routine already has.

What this does not measure #

Not an inconsistent repetition. The third constraint here is consistent: its right-hand side is built to agree with what the first two already say, so the problem has a solution at every ε. A repetition whose right-hand sides disagree makes the problem infeasible, which is a different failure with a different signature — the feasibility residual finally stops being small, and it is the one case where a caller’s usual check works.

Not a rank-deficient B outright. At ε = 0 the constraint block has rank two, the exact reference cannot be computed, and the routines divide by a zero pivot. That case is detectable and detected, and it is the reason the near case matters: a routine that refuses at ε = 0 and returns a confident wrong answer at ε = 10⁻¹² has a discontinuity exactly where the problem does not.

And not a remedy. The obvious one is to factor B with column pivoting, detect the near dependence, and drop or combine the offending row — which changes the problem a caller stated and needs a threshold nothing here supplies. Whether the answer from a truncated constraint set is closer to the intended one than the answer from the full set is not measured, and it is not obvious: the near-dependent constraint does carry a little real information, and dropping it discards that along with the conditioning.

Still open: what the multipliers say #

The one quantity this page has not looked at is the one the saddle-point route produces for free.

A nearly dependent pair of constraints has nearly parallel rows in B, and the multipliers attached to them are the quantities that would blow up — a small change in d requires a large change in λ to absorb it. So ‖λ‖ should rise with κ(B), and it is available from the route that is worst at the problem rather than from the route that is best, which is an awkward place for a diagnostic to live.

Whether it rises proportionally to κ(B), and so carries the same information as the R factor’s diagonal at a much higher price, or whether it says something the conditioning does not, is one sweep. The interesting possibility is the second: the multipliers depend on d as well as on B, so they would distinguish a nearly dependent constraint that happens to be consistent from one that is straining against the others — and the conditioning cannot, because it never looks at a right-hand side.

One line #

Two constraints that nearly repeat each other cost four digits of an answer that satisfies both of them to fifteen, and the quantity that would have said so is the conditioning of a matrix a caller does not think of as one.