Exact arithmetic, and what it costs instead

A prime that divides the answer

A modular elimination reports a singular matrix and is telling the truth — over the field with p elements the matrix is singular. Over the rationals it is not. Nothing in the residue distinguishes the two cases, no quantity is small enough to be suspicious, and the wrong answer is a correct computation of a different question.

Worth reading first: The number that decides nothing · An answer with no error in it · Rank is a decision.

Modular arithmetic buys the exact-arithmetic field its only route with no growth in it: work modulo a prime below 2²⁶, and every quantity from the first step to the last is a machine word. The price is paid in one place, and it is unlike anything else on this site.

Build a six by six integer matrix whose determinant is 30030. Reduce it modulo 3 and eliminate. The elimination reaches a column of zeros and reports that the determinant is zero.

det A = 30030, and its residue at 25 primes — 5 of which divide itEach dot is det A mod p, drawn as a fraction of p so the primes are comparable. A dot on the floor is a prime at which the elimination reports a singular matrix — and the matrix is not singular: its determinant is 30030. 5 of 25 primes here are unlucky, a rate of 20 per cent, which is what it looks like when the primes are small. At the 26-bit primes a modular determinant actually uses, a prime divides a 15-bit determinant with probability of order 2.24·10⁻⁷ — but the consequence is unchanged, and there is no way to see it from one residue. What sees it is the budget: computing more residues than Hadamard's bound requires, and refusing to reconstruct when the extra ones disagree.det mod p, as a fraction of p3713a prime that divides the answerdet A3·10⁴primes swept25unlucky5rate0.2det, in bits15singular mod p is not singularand one residue cannot tell them apart
Fig. 1 The determinant is 30030 = 2·3·5·7·11·13, and its residue at twenty-five small primes. Five dots sit on the floor. Each of them is a correct computation.

That report is true. Over the field with three elements the matrix genuinely is singular: its rows are linearly dependent there, its rank is five rather than six, and every statement the elimination made about it is right. What is wrong is nothing at all inside the computation. What is wrong is that a caller asked about the rationals and was answered about 𝔽₃.

Why this has no floating-point analogue

It is tempting to file this beside the site’s other zero-detection problems and it does not belong there.

Deciding that a zero has arrived is about a quantity that is 10⁻¹⁴ when the theory says zero, and about the judgement required to call it one. Rank is a decision is the same judgement in its most consequential form: a singular value of 10⁻¹² is not zero, calling it zero is a choice, and the choice is where the rank comes from. In both cases there is a small number, somebody has to decide about it, and the essay is about the decision.

Here there is no small number. The residue is the integer 0, exactly, arrived at by a sequence of integer operations none of which rounded. There is nothing to threshold, nothing to be careful about, and no amount of extra precision helps — extra precision is what this route already has infinitely much of. The failure is not in the arithmetic’s accuracy. It is that a different question was answered, correctly.

That makes it the second instance on this site of a class the graph field opened: a computation whose output is right and whose subject is wrong. The first was a partition decided in the last digit, where a correctly computed eigenvector of a degenerate matrix yields a partition that is one of several and the choice among them was made by rounding. Here the substitution is cleaner and more total: not one of several answers, but the answer to a different question, and it does not even have the courtesy to vary with the seed.

How likely it is, and why that is the wrong question

A prime p is unlucky for a matrix when p divides det A. For a determinant D, the number of primes below N that divide it is at most log₂|D|, so among primes near 2²⁶ the proportion that are unlucky is about log₂|D| divided by the number of primes in the range — a probability of order 10⁻⁶ for a forty-bit determinant, and smaller as the primes get larger.

The figure uses small primes deliberately, because at 26 bits the failure would be invisible and this essay would have nothing to draw. Five of the first twenty-five odd primes divide 30030, which is a rate of twenty per cent; at word-sized primes the same matrix would need millions of draws to hit one.

det A = 30030, and its residue at 44 primes — 5 of which divide itEach dot is det A mod p, drawn as a fraction of p so the primes are comparable. A dot on the floor is a prime at which the elimination reports a singular matrix — and the matrix is not singular: its determinant is 30030. 5 of 44 primes here are unlucky, a rate of 11.4 per cent, which is what it looks like when the primes are small. At the 26-bit primes a modular determinant actually uses, a prime divides a 15-bit determinant with probability of order 2.24·10⁻⁷ — but the consequence is unchanged, and there is no way to see it from one residue. What sees it is the budget: computing more residues than Hadamard's bound requires, and refusing to reconstruct when the extra ones disagree.det mod p, as a fraction of p311a prime that divides the answerdet A3·10⁴primes swept44unlucky5rate0.11det, in bits15singular mod p is not singularand one residue cannot tell them apart
Fig. 2 The same sweep taken further. The floor fills slowly and never empties: the six primes dividing the determinant are the only ones that ever land there, and they are the first six.

The temptation is to treat that probability as a defence, and it is not one, for two reasons.

The determinant is not random and neither are the primes. A matrix arriving from a discretisation or a combinatorial construction very often has a determinant with structure in it — a factorial, a power of two, a product of small primes — and the implementations that use fixed prime tables use the same primes every time. A rare event that is deterministic is not rare in the way the arithmetic suggests: it happens to the same inputs every run, for ever.

And the consequence is not a small error. A single unlucky prime among k does not perturb the reconstruction slightly. The Chinese remainder theorem combines residues that are jointly consistent with exactly one value modulo the product, and one wrong residue gives a completely different value — an integer with the right number of digits, none of them right, and no residual to notice it by.

What redundancy buys, and it is not what it looks like

The defence is to compute more residues than the budget requires and compare. It works, and the way it works is worth being precise about, because it is not error correction.

Reconstruct from the first k primes, where k is what Hadamard’s bound required, and reconstruct again from k+1. If no prime was unlucky, both reconstructions give det A and they agree. If one was, the two reconstructions differ, because they are consistent with different values.

So the extra prime does not tell which residue was wrong, and it does not repair anything. It converts an undetected wrong answer into a detected one, which in a field whose entire selling point is that the answer is a certificate is the difference between the route being usable and not. Recovering from it is a retry: drop the primes, take new ones, and the probability that two independent draws are both unlucky is the square of a number that was already small.

10 × 10: three exact determinants agree to the digit, and the float one does notdet A = -77985307738 on a random matrix of order 10, computed three ways with no rounding anywhere: fraction-free elimination, elimination over the rationals, and elimination modulo 2 word-sized primes put back together by the Chinese remainder theorem. The three differ by exactly zero — as integers, not to a tolerance. The floating-point LU of the same matrix returns -7.79853077·10¹⁰, a relative error of 3.913·10⁻¹⁶, which is what a determinant of 37 bits does to a significand of 53. The modular route's entries never exceed 26 bits, which is why the number of primes has to be decided in advance: Hadamard's bound says 43 bits are needed and 2 primes supply them.Bareiss37 bitsrationals36 bitsmodular26 bitsthe answer37 bitsHadamard43 bitswidest number formedrelative error against the exact answerfraction-free (Bareiss)0rational elimination0modular, 2 primes0floating-point LU3.9·10⁻¹⁶det A = -77985307738three routes, one integerand no tolerance anywhere in the comparison
Fig. 3 The reconstruction when nothing goes wrong: three routes, one integer, and the modular one assembled from residues none of which is the answer.

There is a cheaper defence available for a determinant specifically and it is worth naming because it does not generalise. If any single residue is nonzero then det A is nonzero, since a zero determinant has zero residue at every prime. So a nonsingularity certificate needs one lucky prime rather than a full budget — and a singularity claim needs the full reconstruction, because a matrix can be singular modulo every prime anybody tried and not singular. The asymmetry is the useful part: one of the two answers is cheap to prove and the other is not.

The rank moves too, and it does not come back

The determinant is the easy case, because it is one number and a wrong residue makes the whole reconstruction disagree. The rank is worse.

Rank over 𝔽ₚ is computed by the same elimination and it is never larger than the rank over the rationals — a dependency among rows over ℚ is a dependency over 𝔽ₚ, but not conversely. So a modular rank is a lower bound on the true rank, and it is a silent lower bound: the elimination does not report having lost anything, because from its point of view nothing was lost.

Which primes lose rank is not arbitrary either. It is exactly the primes dividing the invariant factors of the matrix — the diagonal of its Smith normal form — which means the failure is a property of the matrix’s structure rather than of the arithmetic. Two essays further on, the invariant factors get computed for their own sake; here they are the answer to which primes are unlucky, and for what.

One 6 × 6 integer matrix, 3 different ranksThe same matrix of integers, its rank computed exactly with no rounding and no threshold anywhere, over the rationals and over the fields with 2, 3, 5, 7, 11, 13, 101, 65537 elements. Over ℚ the rank is 6; over 𝔽2 and 𝔽3 it is 4 and 5. Nothing has been approximated and nothing has been decided by a tolerance — the eliminations are exact in every case, and there is no small number anywhere for a threshold to be compared against. The rank is not a property of the array of numbers; it is a property of the array together with the ring the entries are read in, and a modular algorithm that reports a rank has reported the second one.6𝔽24𝔽35𝔽56𝔽76𝔽116𝔽136𝔽1016𝔽655376rank, by the ring the entries are read inthe rank of one matrixover ℚ6over 𝔽24over 𝔽35over 𝔽56over 𝔽76nothing is rounded hereand the answer still is not a property of the matrix
Fig. 4 One integer matrix, its rank computed exactly over the rationals and over eight finite fields, with no threshold anywhere. Three different answers.

Lifting instead of spreading, and the same trap

There is a second modular strategy that avoids many primes entirely, and it meets the same problem in a form that is easier to get wrong.

Instead of computing modulo p₁, p₂, … and combining, compute modulo a single prime power p, p², p⁴, … by lifting: solve modulo p, use the residual to correct, solve modulo p again, and each round doubles the number of correct p-adic digits. It is the exact-arithmetic analogue of iterative refinement, and the analogy is close enough to be useful and close enough to mislead.

Refinement in floating point converges when the residual is computed accurately enough and the matrix is not too ill conditioned; it is a numerical statement with a numerical failure mode. Lifting converges when the matrix is invertible modulo p, and if it is not, the very first solve fails and nothing after it can recover. There is no gradual degradation and there is no conditioning: the matrix is invertible mod p or the whole route is unavailable.

So the lifting route has to test one prime for luck at the start and can then forget the question, where the spreading route has to worry about every prime it uses. That is a genuine advantage and it comes with a matching disadvantage — one prime is one point of failure, and a prime chosen badly by a fixed table is chosen badly every run.

Which primes are unlucky is a question with an exact answer

The probability argument above is the one usually given and it is the weaker of the two available. There is an exact characterisation, and it is worth stating because it turns a risk into a fact about the matrix.

A prime p is unlucky for the determinant exactly when p divides det A. It is unlucky for the rank exactly when p divides one of the invariant factors — the diagonal entries of the Smith normal form — and the number of ranks lost is the number of invariant factors it divides. The matrix in the hero has invariant factors 1, 1, 1, 1, 1, 30030, so the six primes dividing 30030 each cost it one rank and nothing else does.

That is not an estimate and it does not involve a distribution. It means the set of unlucky primes is a finite, computable, matrix-specific list, and that a routine which knew the invariant factors would know exactly which primes to avoid. Of course computing them is harder than the determinant was — which is the usual shape of these things, and is why the practical defence stays probabilistic.

What it does buy is the right mental model. An unlucky prime is not noise and it is not bad luck in any useful sense; it is one of a handful of specific numbers attached to the matrix, and if the matrix came from a construction with small numbers in it then those numbers are small.

6 × 6: the invariant factors, computed twice by routes with no algorithm in commonThe Smith normal form of a random matrix of order 6 is the diagonal [1, 1, 1, 1, 1, 40098], reached by unimodular row and column operations, with each entry dividing the next — a property the elimination does not give and which has to be enforced and asserted. The same numbers computed from the definition — the gcd of all k × k minors, over the gcd of all (k−1) × (k−1) minors — are [1, 1, 1, 1, 1, 40098], and the two agree entry for entry. Their product is 40098, which is |det A| = 40098. The Hermite form of the same matrix reaches 27 bits against the Smith form's 16, and its transform has determinant -1, so the lattice the rows generate is unchanged.s10 bitss20 bitss30 bitss40 bitss50 bitss615 bitsinvariant factors, in bitstwo routes, and a conserved quantitydet A-4·10⁴Π invariants4·10⁴SNF widest16HNF widest27det of the transform-1an algorithm and a definitionagreeing as integers
Fig. 5 The invariant factors of a random integer matrix: five ones and one large factor carrying the whole determinant. The primes dividing that last factor are the unlucky ones, and there are no others.

What a caller should ask for

The practical upshot is a rule about interfaces, and it applies well outside this field.

A modular routine should not return “the determinant”. It should return the determinant together with the modulus its answer is certain to, because that is what it actually knows — and a caller that wanted a rational answer can then check whether the modulus was large enough while a caller that wanted a residue can use it directly. Returning the number alone throws away the only information that distinguishes a finished computation from an unfinished one.

That is the same complaint the essay about what a method reports makes about an iterative solver returning its own estimated residual as though it were the true one, and the same complaint the certificate essays make about a routine that returns an answer where it could have returned an answer and a proof. A number without its warranty is a number somebody will use outside the range it was true in, which is the failure the essay about a tolerance nobody meets describes from the other end.

The cost of the defence, measured against what it defends

The redundancy is one extra elimination out of k, so it costs 1/k of the route — three per cent on a determinant needing thirty primes, thirty-three on one needing three. That is an unusual cost curve: the check is cheapest exactly where the computation is largest, and most expensive on the small problems where a fraction-free elimination would have been affordable anyway.

Which suggests the rule the field actually follows. Below the size where the modular route wins, fraction-free elimination has no unlucky primes at all and should be used; above it, the redundancy is nearly free and should always be taken. There is no regime in which the modular route without a check is the right answer, and the reason is not caution. It is that the middle regime — where the check is expensive — is exactly the regime where the other route is available.

The refusal

The claim this essay exists to break is the natural one: a nonzero determinant has a nonzero residue at every prime, so a modular elimination reporting zero has found a singular matrix. It sounds like a statement about arithmetic and it is a statement about nothing.

The construction that refuses it is built rather than searched for. Take diag(1, …, 1, 30030) and conjugate it by two unimodular matrices — integer matrices with determinant ±1, so the Smith normal form is untouched and the determinant is exactly 30030. The result is a matrix with no visible structure, entries in the hundreds, nonsingular over ℚ, and singular modulo each of 2, 3, 5, 7, 11 and 13. The assertion is fed that matrix and has to reject it.

Searching random matrices for one would have made the figure a statement about the search. Building it makes the six unlucky primes a stated property of the construction, which is what lets the figure assert that the residue is zero exactly at the primes dividing the determinant rather than merely observing that some are.

The distinction is the same one the site makes about every constructed counterexample. A matrix found by looking is evidence that such matrices exist; a matrix built from a property is evidence about what the property does, and it can carry an assertion that fails if the property stops holding. Here the assertion is an equivalence rather than an example: the residue is zero at p if and only if p divides the determinant, checked at every prime in the sweep, so a change to the construction that broke the correspondence would be caught rather than merely making the picture less striking.

At other settings

det A = 30030, and its residue at 25 primes — 5 of which divide itEach dot is det A mod p, drawn as a fraction of p so the primes are comparable. A dot on the floor is a prime at which the elimination reports a singular matrix — and the matrix is not singular: its determinant is 30030. 5 of 25 primes here are unlucky, a rate of 20 per cent, which is what it looks like when the primes are small. At the 26-bit primes a modular determinant actually uses, a prime divides a 15-bit determinant with probability of order 2.24·10⁻⁷ — but the consequence is unchanged, and there is no way to see it from one residue. What sees it is the budget: computing more residues than Hadamard's bound requires, and refusing to reconstruct when the extra ones disagree.det mod p, as a fraction of p3713a prime that divides the answerdet A3·10⁴primes swept25unlucky5rate0.2det, in bits15singular mod p is not singularand one residue cannot tell them apart
Fig. 6 The same construction at order four, where the determinant is unchanged because it is put there by the diagonal rather than by the size.
det A = 30030, and its residue at 30 primes — 5 of which divide itEach dot is det A mod p, drawn as a fraction of p so the primes are comparable. A dot on the floor is a prime at which the elimination reports a singular matrix — and the matrix is not singular: its determinant is 30030. 5 of 30 primes here are unlucky, a rate of 16.7 per cent, which is what it looks like when the primes are small. At the 26-bit primes a modular determinant actually uses, a prime divides a 15-bit determinant with probability of order 2.24·10⁻⁷ — but the consequence is unchanged, and there is no way to see it from one residue. What sees it is the budget: computing more residues than Hadamard's bound requires, and refusing to reconstruct when the extra ones disagree.det mod p, as a fraction of p3713a prime that divides the answerdet A3·10⁴primes swept30unlucky5rate0.17det, in bits15singular mod p is not singularand one residue cannot tell them apart
Fig. 7 Order eight, thirty primes, and the same five on the floor.
det A = 30030, and its residue at 64 primes — 5 of which divide itEach dot is det A mod p, drawn as a fraction of p so the primes are comparable. A dot on the floor is a prime at which the elimination reports a singular matrix — and the matrix is not singular: its determinant is 30030. 5 of 64 primes here are unlucky, a rate of 7.81 per cent, which is what it looks like when the primes are small. At the 26-bit primes a modular determinant actually uses, a prime divides a 15-bit determinant with probability of order 2.24·10⁻⁷ — but the consequence is unchanged, and there is no way to see it from one residue. What sees it is the budget: computing more residues than Hadamard's bound requires, and refusing to reconstruct when the extra ones disagree.det mod p, as a fraction of p313a prime that divides the answerdet A3·10⁴primes swept64unlucky5rate0.078det, in bits15singular mod p is not singularand one residue cannot tell them apart
Fig. 8 Sixty-four primes swept, of which the unlucky ones are still only the six smallest.
One 8 × 8 integer matrix, 3 different ranksThe same matrix of integers, its rank computed exactly with no rounding and no threshold anywhere, over the rationals and over the fields with 2, 3, 5, 7, 11, 13, 101, 65537 elements. Over ℚ the rank is 8; over 𝔽2 and 𝔽3 it is 6 and 7. Nothing has been approximated and nothing has been decided by a tolerance — the eliminations are exact in every case, and there is no small number anywhere for a threshold to be compared against. The rank is not a property of the array of numbers; it is a property of the array together with the ring the entries are read in, and a modular algorithm that reports a rank has reported the second one.8𝔽26𝔽37𝔽58𝔽78𝔽118𝔽138𝔽1018𝔽655378rank, by the ring the entries are read inthe rank of one matrixover ℚ8over 𝔽26over 𝔽37over 𝔽58over 𝔽78nothing is rounded hereand the answer still is not a property of the matrix
Fig. 9 The rank over eight fields at order eight, where the two small primes lose one and two respectively.
One 4 × 4 integer matrix, 3 different ranksThe same matrix of integers, its rank computed exactly with no rounding and no threshold anywhere, over the rationals and over the fields with 2, 3, 5, 7, 11, 13, 101, 65537 elements. Over ℚ the rank is 4; over 𝔽2 and 𝔽3 it is 2 and 3. Nothing has been approximated and nothing has been decided by a tolerance — the eliminations are exact in every case, and there is no small number anywhere for a threshold to be compared against. The rank is not a property of the array of numbers; it is a property of the array together with the ring the entries are read in, and a modular algorithm that reports a rank has reported the second one.4𝔽22𝔽33𝔽54𝔽74𝔽114𝔽134𝔽1014𝔽655374rank, by the ring the entries are read inthe rank of one matrixover ℚ4over 𝔽22over 𝔽33over 𝔽54over 𝔽74nothing is rounded hereand the answer still is not a property of the matrix
Fig. 10 And at order four, the smallest size the construction exists at.
5 × 5: the invariant factors, computed twice by routes with no algorithm in commonThe Smith normal form of a random matrix of order 5 is the diagonal [1, 1, 1, 1, 22616], reached by unimodular row and column operations, with each entry dividing the next — a property the elimination does not give and which has to be enforced and asserted. The same numbers computed from the definition — the gcd of all k × k minors, over the gcd of all (k−1) × (k−1) minors — are [1, 1, 1, 1, 22616], and the two agree entry for entry. Their product is 22616, which is |det A| = 22616. The Hermite form of the same matrix reaches 24 bits against the Smith form's 15, and its transform has determinant -1, so the lattice the rows generate is unchanged.s10 bitss20 bitss30 bitss40 bitss514 bitsinvariant factors, in bitstwo routes, and a conserved quantitydet A-2.3·10⁴Π invariants2.3·10⁴SNF widest15HNF widest24det of the transform-1an algorithm and a definitionagreeing as integers
Fig. 11 The invariant factors, which are exactly what decides which primes are unlucky.
The determinant's length, and the bound that decides how many primes to use before any of them is chosenSolid: Hadamard's bound on |det A|, the product of the row 2-norms, in bits. Dotted: the determinant's actual length. The bound is what a modular determinant budgets against — at n = 11 the random family needs 47 bits, which is 2 word-sized primes at 26 bits each, and that count is fixed before the first residue is computed. On a Sylvester Hadamard matrix of order 8 the bound is 13 bits and the determinant has 13, so the bound is attained rather than merely satisfied. On the unit triangular family the determinant is 1 at every size while the bound grows, which is the gap the budget pays for: a matrix can have a determinant of one and rows that are long.34567891011110¹10²nbitsthe budget and what it buysrandom, bound47random, actual33primes needed2Hadamard n = 8, bound13Hadamard n = 8, actual13the count is decided by a theorembefore any arithmetic happens
Fig. 12 The budget that says how many residues are needed, and does not say whether any of them is uninformative.
10 × 10: three exact determinants agree to the digit, and the float one does notdet A = -6884259770135383840074327540841 on a wide matrix of order 10, computed three ways with no rounding anywhere: fraction-free elimination, elimination over the rationals, and elimination modulo 5 word-sized primes put back together by the Chinese remainder theorem. The three differ by exactly zero — as integers, not to a tolerance. The floating-point LU of the same matrix returns -6.88425977·10³⁰, a relative error of 9.813·10⁻¹⁶, which is what a determinant of 103 bits does to a significand of 53. The modular route's entries never exceed 26 bits, which is why the number of primes has to be decided in advance: Hadamard's bound says 111 bits are needed and 5 primes supply them.Bareiss103 bitsrationals103 bitsmodular26 bitsthe answer103 bitsHadamard111 bitswidest number formedrelative error against the exact answerfraction-free (Bareiss)0rational elimination0modular, 5 primes0floating-point LU9.8·10⁻¹⁶det A = -6884259770135383840074327540841three routes, one integerand no tolerance anywhere in the comparison
Fig. 13 Four primes on a wide matrix, reconstructing an eighty-bit determinant from four twenty-six-bit residues.
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. 14 The floating-point version of a rank that is not a property of the matrix, which this one is not either, for an unrelated reason.
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.68 wide at rank 8, where the method is at its worst, and 1.004 wide at rank 20, 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.7spread at rank 201best median error0.15widest where the method is worstand the bound does not say so
Fig. 15 A rank as a band rather than a number, which is the shape the float question takes.
Where conjugate gradients certifies that a positive definite matrix is indefiniteEvery matrix in this grid is 16×16 and positive definite by construction — its spectrum is κ^(−i/(n−1)), so the smallest eigenvalue is 1/κ and none of them is negative. A filled cell is a run in which pᵀAp came out non-positive and the iteration produced a direction it would report as a proof of indefiniteness. 18 of the 72 runs did. The region is a staircase whose top edge is at twelve significand bits and whose left edge, at eight, is at κ = 10⁶. Nothing rounded incorrectly anywhere: every comparison was performed exactly as written, on a number that was computed as accurately as the format allows.significand bits10³10⁴10⁵10⁶10⁷10⁸10⁹10¹⁰κ(A)8136460199426314512345105113211154244912624014162024every matrix positive definiteruns producing a false certificate18of runs in total72never above, in significand bits12first κ at eight bits10⁶the comparison was correctand what it proved was not true
Fig. 16 A certificate that is not one, in the field that measures what checking costs.
Four candidate singularity tests on three matrices, at n = 10A table of four quantities on three matrices. 0.1·I at n = 10 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 = 10perfectly 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. 17 Tests for singularity in floating point, none of which have this failure and all of which have another.
Three exact eliminations of the same random matrix, and the widest number each one formsEvery route here returns the same determinant exactly, so there is no error to plot and the y-axis is the length of the intermediates in bits. At n = 11 the answer itself is 33 bits and Hadamard's bound allows 47. Fraction-free elimination never forms a number wider than 35 bits, because each of its intermediates is a minor of the original matrix. Rational elimination that reduces every fraction to lowest terms reaches 35 bits and pays a gcd on every arithmetic operation to stay there. Rational elimination that does not reduce reaches 341287 bits — 9750 times the width of the fraction-free route and 10300 times the width of the answer. All three are correct.34567891011110¹10²10³10⁴10⁵nwidest intermediate, in bitsrationals, not reducedfraction-free · rationals reduced · the answerone answer, three widthsn11answer33Hadamard bound47fraction-free35reduced rationals35unreduced3.4·10⁵the error is zero on every curvethe cost is the length of the numbers
Fig. 18 The growth the modular route removes, which is what it is paying for this risk with.
det A = 30030, and its residue at 25 primes — 5 of which divide itEach dot is det A mod p, drawn as a fraction of p so the primes are comparable. A dot on the floor is a prime at which the elimination reports a singular matrix — and the matrix is not singular: its determinant is 30030. 5 of 25 primes here are unlucky, a rate of 20 per cent, which is what it looks like when the primes are small. At the 26-bit primes a modular determinant actually uses, a prime divides a 15-bit determinant with probability of order 2.24·10⁻⁷ — but the consequence is unchanged, and there is no way to see it from one residue. What sees it is the budget: computing more residues than Hadamard's bound requires, and refusing to reconstruct when the extra ones disagree.det mod p, as a fraction of p3713a prime that divides the answerdet A3·10⁴primes swept25unlucky5rate0.2det, in bits15singular mod p is not singularand one residue cannot tell them apart
Fig. 19 Order five, where the construction puts the same six primes on the floor.
det A = 30030, and its residue at 54 primes — 5 of which divide itEach dot is det A mod p, drawn as a fraction of p so the primes are comparable. A dot on the floor is a prime at which the elimination reports a singular matrix — and the matrix is not singular: its determinant is 30030. 5 of 54 primes here are unlucky, a rate of 9.26 per cent, which is what it looks like when the primes are small. At the 26-bit primes a modular determinant actually uses, a prime divides a 15-bit determinant with probability of order 2.24·10⁻⁷ — but the consequence is unchanged, and there is no way to see it from one residue. What sees it is the budget: computing more residues than Hadamard's bound requires, and refusing to reconstruct when the extra ones disagree.det mod p, as a fraction of p313a prime that divides the answerdet A3·10⁴primes swept54unlucky5rate0.093det, in bits15singular mod p is not singularand one residue cannot tell them apart
Fig. 20 The largest sweep, fifty-four primes, and the same five unlucky ones.
One 5 × 5 integer matrix, 3 different ranksThe same matrix of integers, its rank computed exactly with no rounding and no threshold anywhere, over the rationals and over the fields with 2, 3, 5, 7, 11, 13, 101, 65537 elements. Over ℚ the rank is 5; over 𝔽2 and 𝔽3 it is 3 and 4. Nothing has been approximated and nothing has been decided by a tolerance — the eliminations are exact in every case, and there is no small number anywhere for a threshold to be compared against. The rank is not a property of the array of numbers; it is a property of the array together with the ring the entries are read in, and a modular algorithm that reports a rank has reported the second one.5𝔽23𝔽34𝔽55𝔽75𝔽115𝔽135𝔽1015𝔽655375rank, by the ring the entries are read inthe rank of one matrixover ℚ5over 𝔽23over 𝔽34over 𝔽55over 𝔽75nothing is rounded hereand the answer still is not a property of the matrix
Fig. 21 The rank version of the same substitution, at order five.
One 7 × 7 integer matrix, 3 different ranksThe same matrix of integers, its rank computed exactly with no rounding and no threshold anywhere, over the rationals and over the fields with 2, 3, 5, 7, 11, 13, 101, 65537 elements. Over ℚ the rank is 7; over 𝔽2 and 𝔽3 it is 5 and 6. Nothing has been approximated and nothing has been decided by a tolerance — the eliminations are exact in every case, and there is no small number anywhere for a threshold to be compared against. The rank is not a property of the array of numbers; it is a property of the array together with the ring the entries are read in, and a modular algorithm that reports a rank has reported the second one.7𝔽25𝔽36𝔽57𝔽77𝔽117𝔽137𝔽1017𝔽655377rank, by the ring the entries are read inthe rank of one matrixover ℚ7over 𝔽25over 𝔽36over 𝔽57over 𝔽77nothing is rounded hereand the answer still is not a property of the matrix
Fig. 22 And at order seven, where the two small primes still cost one and two.
4 × 4: the invariant factors, computed twice by routes with no algorithm in commonThe Smith normal form of a unitriangular matrix of order 4 is the diagonal [1, 1, 1, 1], reached by unimodular row and column operations, with each entry dividing the next — a property the elimination does not give and which has to be enforced and asserted. The same numbers computed from the definition — the gcd of all k × k minors, over the gcd of all (k−1) × (k−1) minors — are [1, 1, 1, 1], and the two agree entry for entry. Their product is 1, which is |det A| = 1. The Hermite form of the same matrix reaches 1 bits against the Smith form's 1, and its transform has determinant 1, so the lattice the rows generate is unchanged.s10 bitss20 bitss30 bitss40 bitsinvariant factors, in bitstwo routes, and a conserved quantitydet A1Π invariants1SNF widest1HNF widest1det of the transform1an algorithm and a definitionagreeing as integers
Fig. 23 Every invariant factor 1, so the set of unlucky primes is empty.
The determinant's length, and the bound that decides how many primes to use before any of them is chosenSolid: Hadamard's bound on |det A|, the product of the row 2-norms, in bits. Dotted: the determinant's actual length. The bound is what a modular determinant budgets against — at n = 8 the random family needs 32 bits, which is 2 word-sized primes at 26 bits each, and that count is fixed before the first residue is computed. On a Sylvester Hadamard matrix of order 8 the bound is 13 bits and the determinant has 13, so the bound is attained rather than merely satisfied. On the unit triangular family the determinant is 1 at every size while the bound grows, which is the gap the budget pays for: a matrix can have a determinant of one and rows that are long.345678110¹10²nbitsthe budget and what it buysrandom, bound32random, actual24primes needed2Hadamard n = 8, bound13Hadamard n = 8, actual13the count is decided by a theorembefore any arithmetic happens
Fig. 24 The budget, which counts residues and does not check any of them.
6 × 6: three exact determinants agree to the digit, and the float one does notdet A = -522130948717689337 on a wide matrix of order 6, computed three ways with no rounding anywhere: fraction-free elimination, elimination over the rationals, and elimination modulo 3 word-sized primes put back together by the Chinese remainder theorem. The three differ by exactly zero — as integers, not to a tolerance. The floating-point LU of the same matrix returns -5.22130949·10¹⁷, a relative error of 2.451·10⁻¹⁶, which is what a determinant of 59 bits does to a significand of 53. The modular route's entries never exceed 26 bits, which is why the number of primes has to be decided in advance: Hadamard's bound says 65 bits are needed and 3 primes supply them.Bareiss59 bitsrationals59 bitsmodular26 bitsthe answer59 bitsHadamard65 bitswidest number formedrelative error against the exact answerfraction-free (Bareiss)0rational elimination0modular, 3 primes0floating-point LU2.5·10⁻¹⁶det A = -522130948717689337three routes, one integerand no tolerance anywhere in the comparison
Fig. 25 Three primes on a wide matrix, none of them unlucky.
Three exact eliminations of the same random matrix, and the widest number each one formsEvery route here returns the same determinant exactly, so there is no error to plot and the y-axis is the length of the intermediates in bits. At n = 9 the answer itself is 25 bits and Hadamard's bound allows 37. Fraction-free elimination never forms a number wider than 25 bits, because each of its intermediates is a minor of the original matrix. Rational elimination that reduces every fraction to lowest terms reaches 25 bits and pays a gcd on every arithmetic operation to stay there. Rational elimination that does not reduce reaches 5947 bits — 238 times the width of the fraction-free route and 238 times the width of the answer. All three are correct.3456789110¹10²10³10⁴nwidest intermediate, in bitsrationals, not reducedfraction-free · rationals reduced · the answerone answer, three widthsn9answer25Hadamard bound37fraction-free25reduced rationals25unreduced5947the error is zero on every curvethe cost is the length of the numbers
Fig. 26 What the modular route is buying with this risk.
Recovering 355/113 from its residue: nothing below 2^18, everything aboveThe fraction 355/113 is reduced modulo a prime just above 2^e and then reconstructed from that one residue by the extended Euclidean algorithm. Below e = 18 the reconstruction returns a different fraction or refuses; at and above it, it returns 355/113 exactly, at every e tried. The crossover is not a property of the algorithm: a residue determines a rational uniquely once the modulus exceeds twice the square of the LARGER of the two parts, which here is 2^17.94, and below that there are two fractions with the same residue and both parts small enough. The measured first success is 2^18. There is no band in which the answer is is nearly right — the value returned below the cliff is a perfectly good fraction that is not this one.2^18modulus, as a power of twoa lattice with one short vectornumerator355denominator1132·max(n, d)², bits18first recovered at18moduli tried42below the bound there are two answersand the algorithm cannot prefer one
Fig. 27 The other bound the modular route works to.
The exact solution of a wide system is longer than the systemAx = b with A a wide integer matrix and b chosen so the answer is a genuine vector of fractions. At n = 9 the question occupies 810 bits — 9² entries of 10 bits each — and the answer occupies 1626, a ratio of 2.01. Every entry is a ratio of two determinants by Cramer's rule, checked here by multiplying each one back by det A and comparing it with the determinant of the column-replaced matrix as an integer — so every denominator divides det A and none of them is longer than its 91 bits, while a numerator can be, at 92 here. The answer IS that long, whatever route produced it. That is a floor no exact solver can get under, and it is why the modular route reconstructs the answer rather than eliminating its way to it.345678910²10³nbitsthe answerthe question · det Abits in, bits outn9question810answer1626det A91widest part92the answer is the floorand no route writes it down more cheaply
Fig. 28 The answer a solve has to reconstruct once the residues are in.
det A = 30030, and its residue at 36 primes — 5 of which divide itEach dot is det A mod p, drawn as a fraction of p so the primes are comparable. A dot on the floor is a prime at which the elimination reports a singular matrix — and the matrix is not singular: its determinant is 30030. 5 of 36 primes here are unlucky, a rate of 13.9 per cent, which is what it looks like when the primes are small. At the 26-bit primes a modular determinant actually uses, a prime divides a 15-bit determinant with probability of order 2.24·10⁻⁷ — but the consequence is unchanged, and there is no way to see it from one residue. What sees it is the budget: computing more residues than Hadamard's bound requires, and refusing to reconstruct when the extra ones disagree.det mod p, as a fraction of p311a prime that divides the answerdet A3·10⁴primes swept36unlucky5rate0.14det, in bits15singular mod p is not singularand one residue cannot tell them apart
Fig. 29 Order four, with the same determinant and the same six primes.
det A = 30030, and its residue at 54 primes — 5 of which divide itEach dot is det A mod p, drawn as a fraction of p so the primes are comparable. A dot on the floor is a prime at which the elimination reports a singular matrix — and the matrix is not singular: its determinant is 30030. 5 of 54 primes here are unlucky, a rate of 9.26 per cent, which is what it looks like when the primes are small. At the 26-bit primes a modular determinant actually uses, a prime divides a 15-bit determinant with probability of order 2.24·10⁻⁷ — but the consequence is unchanged, and there is no way to see it from one residue. What sees it is the budget: computing more residues than Hadamard's bound requires, and refusing to reconstruct when the extra ones disagree.det mod p, as a fraction of p313a prime that divides the answerdet A3·10⁴primes swept54unlucky5rate0.093det, in bits15singular mod p is not singularand one residue cannot tell them apart
Fig. 30 Order seven, fifty-four primes swept.
One 4 × 4 integer matrix, 3 different ranksThe same matrix of integers, its rank computed exactly with no rounding and no threshold anywhere, over the rationals and over the fields with 2, 3, 5, 7, 11, 13, 101, 65537 elements. Over ℚ the rank is 4; over 𝔽2 and 𝔽3 it is 2 and 3. Nothing has been approximated and nothing has been decided by a tolerance — the eliminations are exact in every case, and there is no small number anywhere for a threshold to be compared against. The rank is not a property of the array of numbers; it is a property of the array together with the ring the entries are read in, and a modular algorithm that reports a rank has reported the second one.4𝔽22𝔽33𝔽54𝔽74𝔽114𝔽134𝔽1014𝔽655374rank, by the ring the entries are read inthe rank of one matrixover ℚ4over 𝔽22over 𝔽33over 𝔽54over 𝔽74nothing is rounded hereand the answer still is not a property of the matrix
Fig. 31 The smallest size the rank construction exists at.
One 9 × 9 integer matrix, 3 different ranksThe same matrix of integers, its rank computed exactly with no rounding and no threshold anywhere, over the rationals and over the fields with 2, 3, 5, 7, 11, 13, 101, 65537 elements. Over ℚ the rank is 9; over 𝔽2 and 𝔽3 it is 7 and 8. Nothing has been approximated and nothing has been decided by a tolerance — the eliminations are exact in every case, and there is no small number anywhere for a threshold to be compared against. The rank is not a property of the array of numbers; it is a property of the array together with the ring the entries are read in, and a modular algorithm that reports a rank has reported the second one.9𝔽27𝔽38𝔽59𝔽79𝔽119𝔽139𝔽1019𝔽655379rank, by the ring the entries are read inthe rank of one matrixover ℚ9over 𝔽27over 𝔽38over 𝔽59over 𝔽79nothing is rounded hereand the answer still is not a property of the matrix
Fig. 32 And the largest drawn.
4 × 4: the invariant factors, computed twice by routes with no algorithm in commonThe Smith normal form of a wide matrix of order 4 is the diagonal [1, 1, 1, 715250895253], reached by unimodular row and column operations, with each entry dividing the next — a property the elimination does not give and which has to be enforced and asserted. The same numbers computed from the definition — the gcd of all k × k minors, over the gcd of all (k−1) × (k−1) minors — are [1, 1, 1, 715250895253], and the two agree entry for entry. Their product is 715250895253, which is |det A| = 715250895253. The Hermite form of the same matrix reaches 65 bits against the Smith form's 40, and its transform has determinant -1, so the lattice the rows generate is unchanged.s10 bitss20 bitss30 bitss439 bitsinvariant factors, in bitstwo routes, and a conserved quantitydet A-7.2·10¹¹Π invariants7.2·10¹¹SNF widest40HNF widest65det of the transform-1an algorithm and a definitionagreeing as integers
Fig. 33 Invariant factors on wide entries.
12 × 12: three exact determinants agree to the digit, and the float one does notdet A = 647452350646 on a random matrix of order 12, computed three ways with no rounding anywhere: fraction-free elimination, elimination over the rationals, and elimination modulo 3 word-sized primes put back together by the Chinese remainder theorem. The three differ by exactly zero — as integers, not to a tolerance. The floating-point LU of the same matrix returns 6.47452351·10¹¹, a relative error of 3.394·10⁻¹⁵, which is what a determinant of 40 bits does to a significand of 53. The modular route's entries never exceed 26 bits, which is why the number of primes has to be decided in advance: Hadamard's bound says 53 bits are needed and 3 primes supply them.Bareiss40 bitsrationals39 bitsmodular26 bitsthe answer40 bitsHadamard53 bitswidest number formedrelative error against the exact answerfraction-free (Bareiss)0rational elimination0modular, 3 primes0floating-point LU3.4·10⁻¹⁵det A = 647452350646three routes, one integerand no tolerance anywhere in the comparison
Fig. 34 A reconstruction from three residues, none of them unlucky.
Three exact eliminations of the same wide matrix, and the widest number each one formsEvery route here returns the same determinant exactly, so there is no error to plot and the y-axis is the length of the intermediates in bits. At n = 10 the answer itself is 103 bits and Hadamard's bound allows 110. Fraction-free elimination never forms a number wider than 103 bits, because each of its intermediates is a minor of the original matrix. Rational elimination that reduces every fraction to lowest terms reaches 103 bits and pays a gcd on every arithmetic operation to stay there. Rational elimination that does not reduce reaches 754790 bits — 7330 times the width of the fraction-free route and 7330 times the width of the answer. All three are correct.345678910110¹10²10³10⁴10⁵10⁶nwidest intermediate, in bitsrationals, not reducedfraction-free · rationals reduced · the answerone answer, three widthsn10answer103Hadamard bound110fraction-free103reduced rationals103unreduced7.5·10⁵the error is zero on every curvethe cost is the length of the numbers
Fig. 35 The growth this risk is buying freedom from.
The determinant's length, and the bound that decides how many primes to use before any of them is chosenSolid: Hadamard's bound on |det A|, the product of the row 2-norms, in bits. Dotted: the determinant's actual length. The bound is what a modular determinant budgets against — at n = 12 the random family needs 52 bits, which is 3 word-sized primes at 26 bits each, and that count is fixed before the first residue is computed. On a Sylvester Hadamard matrix of order 8 the bound is 13 bits and the determinant has 13, so the bound is attained rather than merely satisfied. On the unit triangular family the determinant is 1 at every size while the bound grows, which is the gap the budget pays for: a matrix can have a determinant of one and rows that are long.3456789101112110¹10²nbitsthe budget and what it buysrandom, bound52random, actual40primes needed3Hadamard n = 8, bound13Hadamard n = 8, actual13the count is decided by a theorembefore any arithmetic happens
Fig. 36 The budget, which counts and does not verify.

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.

Chinese remainder theoremDeterminantExact arithmeticHadamard boundModular arithmeticRank is a decisionUnlucky prime