Two errors, and whose fault they are

The number that decides nothing

The determinant is the first scalar anybody attaches to a matrix and the last one worth consulting. A tenth of the identity has a determinant of 10⁻⁶⁰ and a condition number of exactly one. The Hilbert matrix's determinant stops being right at n = 13 and stops being a number at n = 29, and nothing in between reports either.

Worth reading first: The condition number is an amplifier · An answer that is known.

Ask a room of people who have taken one linear algebra course how to tell whether a matrix can be inverted, and the answer comes back immediately: check whether the determinant is zero. It is the first scalar the subject attaches to a matrix and it arrives with a theorem — A is invertible if and only if det A ≠ 0 — that is correct, complete, and useless.

Useless in a specific sense, which this essay is about. The theorem is a statement about the real numbers, where a quantity is zero or it is not. A floating-point matrix is never exactly singular and never exactly anything else, so the question a computation actually asks is how close to singular, and the determinant answers that question wrongly in both directions, by any margin one cares to name.

Four candidate singularity tests on three matrices, at n = 40A table of four quantities on three matrices. 0.1·I at n = 40 has a condition number of exactly 1 and a determinant of 10⁻⁴⁰; a matrix with κ = 10¹⁰ normalised to |det| = 1 has a determinant of one. The determinant and its nth root give the wrong verdict on both, the smallest singular value gives the wrong verdict on the normalised matrix because it carries the units, and only σ_min/σ_max is right on all three. A test is read as saying *singular* when its number is below 10⁻⁶; green cells agree with the true verdict and red ones do not.does this matrix look nearly singular?green: the test agrees with the truth · red: it does not · the bar under each number is its magnitude, over sixty-two decades|det A||det A|^(1/n)σ_min1/κ = σ_min/σ_max0.1·I at n = 40perfectly conditioned10⁻⁴⁰0.10.11κ = 10¹⁰, |det| = 1nearly singular1110·10⁻⁶10·10⁻¹¹Hilbert at n = 8nearly singular2.7·10⁻³³8.5·10⁻⁵1.1·10⁻¹⁰6.6·10⁻¹¹the two counterexamplesκ of the scaled identity1its determinant10⁻⁴⁰κ of the normalised matrix10¹⁰its determinant1det(cA) = cⁿ det(A)so a determinant carries the units n times over
Fig. 1 Four candidate answers to is this matrix nearly singular, on three matrices chosen so that they disagree. Green means the number gives the right verdict and red means it does not. Drag the size of the scaled identity: only the first column moves, and only its top two rows.

The counterexample is one line long

Take a tenth of the identity, at n = 60.

Its determinant is 10⁻⁶⁰. Its condition number is exactly 1 — not approximately, exactly: every singular value is 0.1, so the ratio of the largest to the smallest is one. Its columns are orthogonal. A solve against it is a division of each component by a tenth, which loses nothing at all, and its inverse is ten times the identity and is exact.

There is no sense in which this matrix is nearly singular. The scalar a reader was taught to consult is 10⁻⁶⁰, which is below every threshold anybody has ever proposed.

The mechanism is not subtle and it is worth naming, because it generalises. det(cA) = cⁿ det(A). The determinant carries the units of the matrix raised to the power of the dimension. Multiply a 60×60 matrix by ten — change nothing about the problem, since the same scaling applied to both sides of Ax = b leaves x alone — and the determinant moves by sixty orders of magnitude.

The previous phase on this site established that the condition number is a choice of units: κ moves under a row scaling that leaves the answer alone. The determinant is worse by a whole dimension. κ is a ratio of two singular values and a uniform scaling divides out of it exactly; det is a product of all n of them and a uniform scaling multiplies it n times over.

Four candidate singularity tests on three matrices, at n = 120A table of four quantities on three matrices. 0.1·I at n = 120 has a condition number of exactly 1 and a determinant of 10⁻¹²⁰; a matrix with κ = 10¹⁰ normalised to |det| = 1 has a determinant of one. The determinant and its nth root give the wrong verdict on both, the smallest singular value gives the wrong verdict on the normalised matrix because it carries the units, and only σ_min/σ_max is right on all three. A test is read as saying *singular* when its number is below 10⁻⁶; green cells agree with the true verdict and red ones do not.does this matrix look nearly singular?green: the test agrees with the truth · red: it does not · the bar under each number is its magnitude, over sixty-two decades|det A||det A|^(1/n)σ_min1/κ = σ_min/σ_max0.1·I at n = 120perfectly conditioned10⁻¹²⁰0.10.11κ = 10¹⁰, |det| = 1nearly singular1110·10⁻⁶10·10⁻¹¹Hilbert at n = 8nearly singular2.7·10⁻³³8.5·10⁻⁵1.1·10⁻¹⁰6.6·10⁻¹¹the two counterexamplesκ of the scaled identity1its determinant10⁻¹²⁰κ of the normalised matrix10¹⁰its determinant1det(cA) = cⁿ det(A)so a determinant carries the units n times over
Fig. 2 The same table at n = 120, where a tenth of the identity has a determinant of 10⁻¹²⁰. Nothing else on the row has moved: the condition number is still one and the matrix is still trivial to solve against.

And the other direction, which is the one that hurts

The converse is worse, because it is the failure mode that lets something through.

Take a matrix with a condition number of 10¹⁰ and divide it by its own geometric mean singular value. That is a uniform scaling, so the condition number is untouched — every singular value moves by the same factor and the ratio does not care. And the determinant is now exactly 1.

A matrix that will destroy ten digits of any answer computed from it, wearing a determinant of one.

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. 3 What the condition number measures and the determinant does not: how much a perturbation of the data is amplified into the answer. Every quantity in this figure is a ratio, which is why a uniform scaling of A does not move it.

The two counterexamples sit at opposite ends of one fact. The determinant is not scale-invariant, and singular is a scale-invariant property. A quantity that moves under a transformation the question does not move under cannot be an answer to the question — which is the previous phase’s finding, arriving in the first scalar the subject introduces.

The obvious repair, and why it is not enough

The scale problem is easy to fix. Take the nth root:

|det A|^(1/n)  =  (σ₁ σ₂ … σₙ)^(1/n)

which is the geometric mean of the singular values, and is invariant under exactly the transformation that broke the determinant. A tenth of the identity has |det|^(1/n) = 0.1 at every size. The scaled identity’s row in the table goes green.

And the repair stops there, because the geometric mean is a mean. One singular value at 10⁻¹⁶ among ninety-nine at 1 moves it by 10^(−16/100), which is 0.69. A matrix whose smallest singular value has vanished — the exact case the question is about — has a scale-free determinant of about seven tenths of its typical singular value, and looks entirely healthy.

A quantity that averages over the spectrum cannot report the smallest member of it, and the smallest member is the question. That sentence is the whole objection and it applies to every determinant-derived quantity there is, because they are all products.

Singular values of a rank-4 matrix with noise of relative size 10⁻⁸Ten singular values on a logarithmic axis. The first four sit near one; the rest sit at the noise level, and the vertical distance between the two groups is the evidence for the rank.1234567891010⁻¹⁸10⁻¹⁵10⁻¹²10⁻⁹10⁻⁶10⁻³1indexsingular valuecutoff, σ₁ · 10⁻¹⁰numerical rank 10gap 8.2·10⁶an opiniontrue rank 410×10, built with 4 nonzero valuesrank is a decision
Fig. 4 The singular value spectrum, which is what a rank question is actually about: the answer is a gap, and a gap is invisible to any product of all the values.
A randomised rank-k solve, 4 seeds a rankRelative error against the rank kept, on a logarithmic vertical axis, with a vertical bar at each rank spanning the seeds. The band is 1.84 wide at rank 8, where the method is at its worst, and 1.022 wide at rank 24, where it is at its best. The same computation on the same data returns a different answer each time.0481216202428321rank keptrelative errormedianthe answer movesspread at rank 81.8spread at rank 241best median error0.14widest where the method is worstand the bound does not say so
Fig. 5 And the band inside which the rank is a decision rather than a fact. Nothing that averages over the spectrum can locate this band, because the band is defined by one singular value.

The number that does answer the question is σ_min/σ_max — scale-free, and a statement about the worst direction rather than the typical one. It costs an SVD, or an estimate of one, and an estimate that can be fooled is about what that estimate costs in honesty.

The determinant is also not a number

Everything above is about what the determinant means. There is a separate and equally serious problem, which is that at ordinary sizes it cannot be represented at all.

A determinant is a sum of n! products of n entries. It therefore scales like the nth power of the entries, and n is not small. A 200×200 matrix of entries around 10 has a determinant of about 10²⁰⁰ times whatever the structure contributes, and binary64 stops at 1.8·10³⁰⁸.

The Hilbert matrix goes the other way, and it goes there fast. Its determinant is known in closed form — the entries are 1/(i+j−1), which is a Cauchy matrix, and Cauchy determinants have a formula — so this site can compute the exact answer in BigInt rationals and compare.

log₁₀|det H_n| by three routes, to n = 30The exact rational determinant of the Hilbert matrix falls to 10^-518 at n = 30. The accumulated logarithm of the pivots follows it to n = 13 and then departs, reaching 10^-352 — wrong by 166 orders of magnitude and still an ordinary-looking number. The product of the pivots is a separate curve that ends at n = 29, where it underflows to exactly zero.26101418222630-550-450-350-250-150-50nlog₁₀|det H_n|exact rationalΣ log|u_ii|Π u_iithree routes, one theoremexact at the largest n-518accumulated logarithm-352decades of disagreement166smallest pivot at that n1.2·10⁻¹⁷every pivot is a normal numberat every size on this axis
Fig. 6 log₁₀|det H_n| by three routes. The exact rational determinant, the accumulated logarithm of the pivots, and the product of the pivots — which ends, because it underflows to exactly zero.

At n = 28 the pivot product is 1.04·10⁻³²¹, which is a subnormal: it is inside the format only because binary64 degrades gracefully below 2.2·10⁻³⁰⁸, and it has already lost most of its significand. The numbers below the smallest one is about what that costs.

At n = 29 it is exactly zero.

And every pivot at n = 29 is a perfectly ordinary number. The smallest is 1.2·10⁻¹⁷, the largest is 1, the factorisation completed, no assertion is false and nothing in the computation went wrong. What left the format is the answer, not the arithmetic.

The exact solution of a 13×13 Hilbert system beside the computed oneTwo columns of numbers: the exact answer, which is the integers one to thirteen, and the answer double-precision elimination returns, with the number of correct digits beside each.H13 x = b, b formed exactly so that x = (1, 2, …, 13)12345678910111213exact1.00002.00003.00133.97865.18864.993010.46720.048921.2657-2.575319.21478.906013.5113computed6.6 correct digits4.8 correct digits3.4 correct digits2.3 correct digits1.4 correct digit0.8 correct digitno correct digitsno correct digitsno correct digitsno correct digitsno correct digits0.6 correct digit1.4 correct digitbackward error2.2·10⁻¹⁷κ = 1.7·10¹⁸The algorithm solved a neighbouring problem perfectly. That problem's answer is this one.right-hand side built in BigInt rationalsthe truth is known
Fig. 7 The Hilbert system whose exact answer is known in closed form. The same closed form gives the exact determinant, which is why the sweep above has a third curve rather than two curves and a hope.

So accumulate the logarithm, which repairs the wrong thing

The fix is standard and it is the shape of this whole phase: do not form the object. The pivots are in hand; take their logarithms and add.

log|det A|  =  Σ log|u_ii|,   sign = parity of the interchanges × signs of the pivots

One logarithm per pivot on top of a factorisation that was going to happen anyway, and the result cannot overflow until the logarithm does, which is never. At n = 30 it returns −1,194, and −1,194 is an entirely ordinary double.

Except that at n = 30 the true value is −1,194 and the computed value is −811.

log₁₀|det H_n| by three routes, to n = 22The exact rational determinant of the Hilbert matrix falls to 10^-274 at n = 22. The accumulated logarithm of the pivots follows it to n = 13 and then departs, reaching 10^-225 — wrong by 49 orders of magnitude and still an ordinary-looking number. The product of the pivots is a separate curve that ends at n = 22, where it underflows to exactly zero.25811141720-300-250-200-150-100-500nlog₁₀|det H_n|exact rationalΣ log|u_ii|Π u_iithree routes, one theoremexact at the largest n-274accumulated logarithm-225decades of disagreement49smallest pivot at that n4·10⁻¹⁷every pivot is a normal numberat every size on this axis
Fig. 8 The same sweep stopped at n = 22, where the two float curves have already separated from the exact one by more than a hundred orders of magnitude in the determinant while still looking like two smooth curves.

The accumulated logarithm departs from the truth at n = 13, and the product underflows at n = 29. Sixteen sizes apart, and the accuracy limit comes first. At n = 20 the computed determinant is 10³⁰·⁶ times the true one — a perfectly ordinary-looking negative number, with nothing on the page to say it is wrong by thirty orders of magnitude.

This is the finding this essay was written to reach and it was not the one expected. The substitution repairs representability and does nothing for conditioning, and the conditioning gives out first by a wide margin. A reader who takes “accumulate the logarithm” as the fix has bought sixteen sizes of an answer that is absent and inherited sixteen sizes of an answer that is wrong.

And the two float routes agree with each other

The trap has a second floor. The site’s standing habit is two routes to every number, and the natural second route to a log-determinant is Σ log σᵢ from the singular values — a completely different algorithm, sharing no arithmetic with the pivots.

At n = 20 the pivots give −448.55 and the singular values give −449.09. They agree to three digits. The exact answer is −518.95.

Two independent computations, agreeing with each other and not with the truth, because both are computing the determinant of the same rounded matrix and the rounded matrix is the problem. The Hilbert matrix’s entries are 1/3, 1/5, 1/7 — none of them representable — so what the routine is handed at n = 20 is a matrix whose determinant genuinely is about 10⁻⁴⁴⁹, and it computes it correctly.

That is the sharpest version of this site’s own sentence. The algorithm is blameless. The problem is sensitive by 10¹⁸. And no second route can find it, because the two routes share the input.

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: κ · u30×30, 20 seeds per κ; dashed is the worstthe problem worsens, not the method
Fig. 9 The two errors a computation can and cannot see. The determinant sweep is this picture with the backward error at rounding and the forward error at thirty orders of magnitude, on a quantity nobody prints a residual for.

The one exact part, which is also wrong

Something in a determinant does not round. The sign is a parity of row interchanges multiplied by the signs of the pivots: comparisons and an integer count, with no arithmetic in them at all.

The Hilbert matrix is positive definite. Every pivot of its exact unpivoted elimination — which is Cholesky’s pivot sequence — is positive, and det H_n > 0 at every n. That is a theorem.

In binary64 the fourteenth pivot is negative, and the reported sign of det H₁₄ is −1. It is wrong again at 15, 18, 19, 24, 25 and 26.

No rounding happened in the sign. What happened is that an exact function was applied to an inexact input, and it returned a wrong answer with no error term attached to it. An inaccurate number carries its own warning, since one can ask how inaccurate. A wrong sign is a clean, confident, integer falsehood.

The growth factor of a 12×12 elimination, against the condition number of the matrixThree curves against κ. Cholesky's growth factor is exactly 1 at every condition number drawn — the elimination never produces an entry larger than the matrix already had. Partial pivoting on the same matrices reaches the same growth, and reaches it by making up to 9 row interchanges, each of which destroys the symmetry that was the reason to use a symmetric factorisation. The bound the general theory allows, 2^11 = 2048, is drawn above them both.10¹10³10⁵10⁷10⁹10¹¹110¹10²10³10⁴condition number of the matrixgrowth factorbound 2^11partial pivotingCholeskyno pivot to gain fromCholesky growth, every κ1Cholesky interchanges0partial pivoting, worst9the bound, 2^112048both eliminations reach the same growthand only one of them had to swap to get there
Fig. 10 The symmetric elimination whose pivots are the ones in question. Its growth factor is exactly one, which is why nothing in the factorisation looks unusual at the size where its pivots start coming out negative.
How often Cholesky still calls a 12×12 matrix positive definite, against λ_min/λ_max in units of the format's own roundoffThree curves, one per precision, of the share of 24 seeded matrices on which the factorisation succeeds. Measured in units of each format's unit roundoff the three lie almost on top of one another, with the edge — the smallest ratio at which every seed succeeds — at 1.0, 1.0, 1.8 times u. The absolute thresholds are 6·10⁻⁸, 2.3·10⁻¹⁰ and 2·10⁻¹⁶: nine orders of magnitude apart, and the same number in the format's own units.10⁻³10⁻²10⁻¹110¹10²10³10⁴00.250.50.751λ_min / λ_max, in units of the format's own ushare of seeds that succeed24 bits32 bits53 bitsone threshold, three formats24-bit edge, absolute6·10⁻⁸32-bit edge, absolute2.3·10⁻¹⁰53-bit edge, absolute2·10⁻¹⁶in units of u, at 53 bits1.8a yes-or-no question with a precision in itand a coin flip three decades below the edge
Fig. 11 And the edge at which a definiteness test stops being a fact about the matrix. The same threshold governs the sign above: below it, whether a pivot comes out positive is a property of how the matrix was rounded when it was stored.

What the determinant is exactly right for

This is not an argument that the determinant is a bad idea, and two of its jobs are untouched.

Volume. |det A| is the factor by which A multiplies volume, and that is not an approximation to anything — it is what the determinant is. Hadamard’s inequality says the volume of a parallelepiped is at most the product of its edge lengths,

|det A|  ≤  ‖a₁‖ ‖a₂‖ … ‖aₙ‖

with equality exactly when the edges are orthogonal. Written as a ratio it is scale-free — both sides carry the same units — and it is a genuinely useful number: a volume efficiency between 0 and 1 that is 1 for an orthogonal matrix and falls as the columns line up. On random 8×8 Gaussian matrices it averages 0.037.

Orientation. The sign, where the input is exact. A permutation matrix, an integer matrix, a rotation given by its entries: for those the sign is a fact and the whole objection above evaporates, because the objection was about an exact function of an inexact input.

Loss of orthogonality against condition number, in binary32A log–log plot of the norm of Q-transpose-Q minus the identity against condition number. Classical Gram–Schmidt rises steeply, modified Gram–Schmidt rises gently, and Householder is flat.110²10⁴10⁶10⁸10¹⁰10¹²10⁻¹⁷10⁻¹⁴10⁻¹¹10⁻⁸10⁻⁵10⁻²10¹condition number κ(A)‖QᵀQ − I‖classicalmodifiedHouseholderκ²uκu8×8, eight seeds per κ, binary32all three reconstruct A
Fig. 12 Orthogonality measured rather than asserted, which is what the Hadamard ratio is a determinant-flavoured version of: 1 exactly when the columns are orthogonal, and falling as they stop being.
QᵀQ from classical Gram–Schmidt and from Householder on the 8×8 Hilbert matrixTwo eight-by-eight tables of QᵀQ. The upper one has ones on the diagonal and entries as large as one off it; the lower one is the identity to three decimal places everywhere.A = H8 · κ = 1.5·10¹⁰ · both factorisations reconstruct A to 6·10⁻¹⁷the diagonal is 1 in both — every column is a unit vector either way1.000000000001.000000000001.000000000001.000000000001.00000.002-0.002000001.0000.125-0.13300000.0020.1251.000-1.0000000-0.002-0.133-1.0001.000classical Gram–Schmidt1.000000000001.000000000001.000000000001.000000000001.000000000001.000000000001.000000000001.000Householderclassical ‖QᵀQ − I‖1.4Householder ‖QᵀQ − I‖1.4·10⁻¹⁵largest off-diagonal 1 against 3.1·10⁻¹⁶length is not angle
Fig. 13 And the same statement entry by entry. A determinant compresses all of this into one number, which is the operation this essay is about the cost of.

Where it is used anyway, and what to use instead

Four places where a determinant appears in code, and what belongs there:

A singularity test. Use σ_min/σ_max, or an estimate of it — every library’s rcond — and never abs(det(A)) < tol. LAPACK does not compute a determinant at all; xGECON estimates the reciprocal condition number and that is the number a driver checks.

A multivariate normal density. The log-likelihood needs log|det Σ|, not det Σ, and Σ is positive definite, so the honest route is 2 Σ log L_ii from a Cholesky factor. Every statistics library that computes a Gaussian density this way is doing exactly the substitution above — and inherits exactly the accuracy limit, on covariance matrices that are routinely ill-conditioned.

A characteristic polynomial. det(λI − A) as a route to eigenvalues is the classic bad idea and the site has drawn what happens: the polynomial’s coefficients are catastrophically worse conditioned than the eigenvalues themselves, and the roots of a computed polynomial are not the eigenvalues of the matrix.

A change of variables. Here the determinant is the answer, it is usually 2×2 or 3×3, the entries are exact, and none of this applies.

How far a perturbation of size ε moves an eigenvalue, 8×8A log–log plot of eigenvalue movement against perturbation size. The symmetric case lies on a line of slope one; the non-symmetric case lies on a line of slope one eighth, and at a perturbation of ten to the minus sixteen it has already moved by a hundredth.10⁻¹⁶10⁻¹³10⁻¹⁰10⁻⁷10⁻⁴10⁻¹10⁻¹⁶10⁻¹³10⁻¹⁰10⁻⁷10⁻⁴10⁻¹size of the perturbation ‖δA‖how far the eigenvalues moveJordan block, ε^(1/8)symmetric, ≤ ‖δA‖rounding error alone moves it to 10⁻²six seeds per symmetric point; Jordan is closed formsymmetry beats precision
Fig. 14 The eigenvalues, which are what a characteristic polynomial is usually a route to. Going through the polynomial’s coefficients puts a badly conditioned intermediate between the matrix and the answer, which is this essay’s subject in a different field.
Backward error, forward error and the condition number, with measured valuesTwo boxes at the top — the problem posed and the nearby problem the algorithm answered exactly — and two answers below them, with the distances between all four labelled by numbers from a Hilbert solve.the problem you posedA = H10b = A·(1, 2, …, 10)the problem it answered exactlyA + δA, b + δb‖δ‖ / ‖A‖ = 2.3·10⁻¹⁷the answer you wantedx = (1, 2, …, 10), exactlythe answer you gotx̂, wrong by 2.7·10⁻⁴ relativebackward error 2.3·10⁻¹⁷forward error 2.7·10⁻⁴κ = 1.6·10¹³κ · η = 3.6·10⁻⁴, and the measured forward error is 2.7·10⁻⁴.The algorithm is not at fault. The problem is.H10, LU with partial pivotingresidual and error differ
Fig. 15 The site’s identity, with the determinant nowhere in it. Forward error, backward error, and the condition number that multiplies them — none of the three is a product over the spectrum, and that is not a coincidence.

What is worth carrying

The determinant answers a question about volume and is asked a question about conditioning. Those are different questions, and the second one is scale-free while the determinant is not, by a whole power of the dimension.

Making it scale-free does not make it right, because it is still a mean. The geometric mean of the singular values cannot see one small singular value among ninety-nine ordinary ones, and one small singular value among ninety-nine ordinary ones is the entire subject.

And two limits arrive, sixteen sizes apart, with the wrong one first. The product underflows at n = 29 and the answer stopped being right at n = 13. The repair everybody knows — accumulate the logarithm — buys the second limit and nothing at all for the first, and the number in between is wrong by thirty orders of magnitude while looking entirely ordinary.

The next essay takes the same subject into the one place a determinant is still taught as a method: a rule that is correct and unusable.

What links here

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

Reads more easily once this is understood

Essays that name this one as worth reading first.

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.

Condition numberDeterminantExact ground truthHilbert matrixLU factorisationOrthogonalityRankSingular valuesSubnormalUnderflow