Number Sieves · Volume 9
Number Sieves — Volume 9 — What a Sieve Can Solve
Worked problems from easy to hard — primes, Sun Zi's remainders, Fermat's squares, Lehmer's quadratic forms, Mersenne numbers and pseudosquares — with the moduli, the survivors and the running time on three machines for each, and the cases where sieving fails
9.1 About This Volume
The earlier volumes describe sieves as machines: loops of bicycle chain with pins in some links, gears drilled with holes that let light through, rings of bits in a microcontroller. This volume is about the work those machines were built to do. It is a ladder of problems, from one a child can do on squared paper to ones that kept D. H. Lehmer’s machines busy for hours, and each rung is worked through in full.
Every problem is set up the same way, because every sieve problem is the same problem in disguise. There is a long row of candidates, numbered 0, 1, 2, 3, … . There is a set of moduli — the chain lengths or gear sizes. For each modulus there is a list of allowed positions, worked out in advance by hand (or by a host computer). A candidate survives if, and only if, it sits on an allowed position for every modulus at once. The machine’s job is to step through the candidates and stop at the survivors. The mathematician’s job — the part that takes thought — is to turn a question about primes, factors or equations into that list of allowed positions.
For each problem the volume gives:
- how the question becomes a sieve (what the candidate is, and what “allowed” means for each modulus);
- how many moduli are needed, and how many survivors to expect;
- a small table of residues or survivors, computed for this volume in Python and checked;
- how long the run takes on three machines: the 1928 bicycle-chain sieve at about 50 candidates per second, the 1932 photoelectric gear sieve at about 5,000 per second, and kbsriram’s 2026 ATtiny412 array at about 500,000 per second. Those three rates are the ones kbsriram uses to compare the machines in his Hackaday.io project log, where he reran Lehmer’s published problems and matched the published times to within half a percent.
Times are simple arithmetic — candidates stepped divided by rate — and they ignore the time spent setting up the moduli, which on the chain machine meant re-pinning nineteen chains by hand. The worked numbers in this volume were all computed and checked in Python; the scripts that draw the figures are in 03-outputs/diagrams/diagrams_c.py. The last section is about failure: the problems a sieve cannot do, and why.
9.2 Rung One: All the Primes Below 100
The first problem is the one the series started with. Volume 2 works the Sieve of Eratosthenes on the numbers 1 to 100 in full, and there is no need to repeat it; what matters here is to see it as the simplest possible case of the machine problem.
The candidates are 2 to 100. There are four moduli: 2, 3, 5 and 7, the primes up to √100 = 10. For each modulus, every position is allowed except the multiples — and the prime itself is spared. Run in that order, the count of survivors falls like this:
Table 1 — The candidates are 2 to 100. There are four moduli: 2, 3, 5 and 7, the primes up to √100 = 10. For each modulus, every position is allowed except the multiples — and the prime itself is spared. Run in that order, the count of survivors falls like this
| Modulus applied | Removed | Survivors left |
|---|---|---|
| (start) | — | 99 |
| 2 | 49 | 50 |
| 3 | 16 | 34 |
| 5 | 6 | 28 |
| 7 | 3 | 25 |
The 25 survivors are the primes below 100. Notice two things that will matter later. First, each modulus does less work than the one before: 2 removes half of what remains, 7 removes only three numbers (49, 77 and 91). Second, the sieve here is exclusive — a candidate survives by avoiding one residue class (zero) for each modulus. Lehmer’s machines worked the other way round: they were inclusive sieves, in which each modulus lists the residues a candidate is allowed to have, usually about half of them. The arithmetic is the same; the bookkeeping is inverted.
Machine time. 99 candidates: 2 seconds on the chain sieve, 20 milliseconds on the gear sieve, 0.2 milliseconds on the ATtiny array. Nobody ever needed a machine for this rung — which is exactly why it is the right place to start.
9.3 Rung Two: Sun Zi’s Remainder Puzzle
The oldest surviving problem that is unmistakably a sieve problem comes from the Sunzi Suanjing (孫子算經, “Master Sun’s Mathematical Manual”), a Chinese arithmetic textbook usually dated to between the third and fifth centuries CE. As translated in the Wikipedia article on the Chinese remainder theorem, it reads: “There are certain things whose number is unknown. If we count them by threes, we have two left over; by fives, we have three left over; and by sevens, two are left over. How many things are there?”
As a sieve. The candidates are the whole numbers 0, 1, 2, … . Three moduli: 3, 5 and 7. Each allows exactly one position:
Table 2 — As a sieve. The candidates are the whole numbers 0, 1, 2, … . Three moduli: 3, 5 and 7. Each allows exactly one position
| Modulus | Allowed residue | Allowed candidates, first few |
|---|---|---|
| 3 | 2 | 2, 5, 8, 11, 14, 17, 20, 23, 26, … |
| 5 | 3 | 3, 8, 13, 18, 23, 28, 33, … |
| 7 | 2 | 2, 9, 16, 23, 30, 37, 44, … |
Because 3, 5 and 7 have no common factor, the pattern of all three repeats every 3 × 5 × 7 = 105 candidates, and within each block of 105 exactly one number passes all three. Stepping through 0 to 104: 35 numbers pass the 3-ring; of those, seven pass the 5-ring (8, 23, 38, 53, 68, 83 and 98); of those, one passes the 7-ring — 23. Every other answer is 23 plus a multiple of 105.
This is the “three-loop sieve” in its purest form. On a chain machine it would be three loops of chain, 3, 5 and 7 links long, each with a single pin; the machine would stop at the first step where all three pins line up under their contacts. Sun Zi’s own method, and the general method that Qin Jiushao published in 1247, solve the problem by calculation rather than search. A sieve does not care: it simply steps until the pins align.
How many moduli, how many survivors. Three moduli, each passing one residue in its own length: the fraction of candidates that survive is 1/3 × 1/5 × 1/7 = 1/105, so one survivor per cycle, and that survivor is the answer. The general lesson: when the moduli are coprime, the survivors of an inclusive sieve are spaced exactly according to the product of the moduli, and a sieve with k allowed residues per modulus finds them all in one pass.
Machine time. At most 105 candidates (to the first answer, 24): half a second on the chain machine and well under a millisecond on the ATtiny array.
9.4 Rung Three: Fermat’s Difference of Squares
The step from remainder puzzles to factoring is Fermat’s. In 1643 Pierre de Fermat explained, in letters to Marin Mersenne, that an odd number which is not a square can be written as a difference of two squares in as many ways as it has factor pairs — and so a factor can be found by finding the squares. His worked example was 2,027,651,281, which he split as 44,021 × 46,061 (Daniel Perrin’s notes for a 2018 Bibliothèque nationale de France lecture; the number is also recorded on the Prime Curios page for 2027651281). The idea is that if
N = a² − b² = (a − b)(a + b),
then a − b and a + b are factors. Start with a just above √N and step upward, testing whether a² − N is a perfect square. For Fermat’s own number √N is just over 45,029, and the square turns up at a = 45,041, where a² − N = 1,020²: twelve steps, checked in Python. Volume 4 works the method by hand for N = 7,663. Here it is applied to a new and larger number, with the sieve doing the testing.
The number. N = 10,012,997. Its square root is 3,164.3, so the candidates start at a = 3,165.
As a sieve. The candidate is a. The test “a² − N is a perfect square” is too expensive to run at every step, but a necessary condition is cheap: for every modulus m, a² − N must be a square modulo m. For each m, the allowed residues of a can be listed once, in advance. Here are the lists for the first seven moduli, as computed in Python:
Table 3 — As a sieve. The candidate is a. The test "a² − N is a perfect square" is too expensive to run at every step, but a necessary condition is cheap: for every modulus m, a² − N must be a square modulo m. For each m, the allowed residues of a can be listed once, in advance. Here are the lists for the first seven moduli, as computed in Python
| Modulus m | N mod m | Allowed residues of a | Fraction allowed |
|---|---|---|---|
| 16 | 5 | 3, 5, 11, 13 | 4/16 |
| 9 | 2 | 0, 3, 6 | 3/9 |
| 5 | 2 | 1, 4 | 2/5 |
| 7 | 1 | 1, 3, 4, 6 | 4/7 |
| 11 | 5 | 3, 4, 5, 6, 7, 8 | 6/11 |
| 13 | 7 | 2, 4, 6, 7, 9, 11 | 6/13 |
| 17 | 14 | 1, 4, 7, 8, 9, 10, 13, 16 | 8/17 |
Prime powers like 16 and 9 are used instead of 2 and 3 because squares modulo a prime power carry more information than squares modulo the prime alone — which is why Lehmer’s machines used prime-power moduli such as 9, 25, 49 and 121 in place of 3, 5, 7 and 11 (Lehmer 1934). Which powers of 2 and 3 the gear sieve carried is not settled; Volume 7 sets out the evidence.
Running it. Suppose the factors are known to lie within a factor of 2.5 of each other, which bounds the search: the answer must come by a = 3,501. That gives 337 candidates, a = 3,165 to 3,501. Applying the moduli one at a time:
Table 4 — Running it. Suppose the factors are known to lie within a factor of 2.5 of each other, which bounds the search: the answer must come by a = 3,501. That gives 337 candidates, a = 3,165 to 3,501. Applying the moduli one at a time
| After modulus | Survivors |
|---|---|
| (none) | 337 |
| 16 | 85 |
| 9 | 29 |
| 5 | 12 |
| 7 | 7 — a = 3,189, 3,219, 3,261, 3,291, 3,429, 3,459, 3,501 |
| 11 | 5 — 3,219, 3,261, 3,429, 3,459, 3,501 |
| 13 | 2 — 3,261, 3,501 |
| 17 | 1 — 3,501 |
At a = 3,501, a² − N = 12,257,001 − 10,012,997 = 2,244,004 = 1,498². So
N = (3,501 − 1,498)(3,501 + 1,498) = 2,003 × 4,999.
Seven moduli reduced 337 candidates to one, and the one was the answer. The order in which moduli are listed makes no difference to the final survivors — a machine applies them all at once — but it shows that each modulus removes roughly half of what remains, which is the rule of thumb for this kind of sieve: an odd prime p passes (p + 1)/2 of its p residues when N is a square modulo p, and (p − 1)/2 when it is not.
Machine time. 337 candidates: about 7 seconds on the chain sieve, 67 milliseconds on the gear sieve, under a millisecond on the ATtiny array. The same number will turn up again in Volume 10, factored a completely different way by the quadratic sieve.
9.5 Rung Four: Quadratic Forms, the Way Lehmer Posed Problems
Fermat’s method works well only when the two factors are close together. Lehmer’s machines were rarely used that way. The problems he and his father published were posed in terms of quadratic forms — expressions like x² + y², x² + 7y² or x² − Dy² — and the machine’s task was to find the values of y for which N − Dy² (or N + Dy²) is a perfect square. The point is a theorem going back to Euler: a number that can be written in the form x² + Dy² in two different ways is composite, and the two representations give away a factor. A number that can be written that way in exactly one way (with the right conditions) is prime. So a sieve that finds every representation both factors numbers and proves them prime.
9.5.1 A small case: two sums of two squares
Take n = 1,022,117. Searching all a up to √n for which n − a² is a perfect square finds exactly two representations:
1,022,117 = 271² + 974² = 314² + 961².
Euler’s rule says that a common factor is hidden in the cross-products. Taking (a, b) = (271, 974) and (c, d) = (314, 961):
gcd(n, ad + bc) = gcd(1,022,117, 566,267) = 1,013, and gcd(n, ad − bc) = 1,009.
So 1,022,117 = 1,009 × 1,013. Both factors are 1 more than a multiple of 4, which is why the number has two representations at all; a prime of the form 4k + 3 can never be a sum of two squares.
9.5.2 A 1928 problem: is 4,500,383,424,373 prime?
D. H. Lehmer’s 1928 American Mathematical Monthly paper on the bicycle-chain machine, “The mechanical combination of linear forms,” includes a problem of exactly this kind. To finish factoring the number 3,173,868,090,153,260² + 1, the Lehmers needed to know whether the factor N = 4,500,383,424,373 was prime. Their method: search every representation N = a² + b². If there is exactly one, N is prime.
As kbsriram reconstructs the run from the paper in his problem driver, a congruence argument fixed a to two progressions, so the machine was set up twice. In the second run, a = 5x + 3, and the candidate is x. For each of the nineteen chains (64, 27, 25, 49, 22, 26 and the primes 17 to 67), the allowed links are the x for which N − (5x + 3)² is a square modulo the chain length. The machine stopped at x = 157,044, which gives a = 785,223 and
N = 785,223² + 1,970,738²,
checked in Python. The first run found nothing. One representation means N is prime, and with it the factorisation of 3,173,868,090,153,260² + 1 is complete: in Python it comes out as 17² × 1,721 × 4,500,383,424,373².
How many moduli, how many survivors. Nineteen chains, but two of them (64 and 25) did almost no filtering because the step size already absorbed those factors. kbsriram’s rerun at the same nineteen moduli found one false stop before the answer, at x = 102,151 (a = 510,758): there N − a² is a quadratic residue of every prime up to 67 and is still not a square. Its smallest non-residue prime is 73 — one chain longer than the 1928 machine had. The paper says the machine “stopped only once”; the reconstruction finds two stops, which kbsriram flags as a departure from the published account.
Machine time. 157,044 candidates to the answer: 52 minutes on the chain sieve (kbsriram’s board, stepping at the same 50 per second, took 52 minutes 37 seconds), 31 seconds on the gear sieve, a third of a second on the ATtiny array.
9.5.3 A 1933 problem: x² + 7y²
In February 1933 Lehmer announced five new factorisations made in “a month’s experiment” with the photoelectric gear sieve, during which “the machine ran only a few hours altogether” (D. H. Lehmer, “Some new factorizations of 2ⁿ ± 1,” Bulletin of the AMS 39 (1933) 105–108). One of them shows how he used forms to prove primality. The large factor of 2⁹⁵ + 1 is 3,011,347,479,614,249,131. To prove it prime, Lehmer used a theorem that needs a large prime factor of N − 1. N − 1 = 570 × 5,283,065,753,709,209, and that 16-digit cofactor M had to be factored first.
Because −7 is a quadratic residue of M, M can be written as x² + 7y². The gear machine searched for representations and found two:
M = 66,855,539² + 7 × 10,779,628² = 40,923,451² + 7 × 22,704,112²,
as printed in the paper and checked in Python. The cross-product rule then gives
gcd(M, 66,855,539 × 22,704,112 + 40,923,451 × 10,779,628) = 59,957,
so M = 59,957 × 88,114,244,437, both prime. kbsriram’s driver for this problem (prob_m16.py) sets y = 4k (forced because M ≡ 1 mod 8), so the candidate k runs from 0 to 6,868,058, and computes that the thirty gears of the 1932 machine together pass about 1 candidate in 100 million — so over the whole range, fewer than 0.07 false stops are expected. The two real representations are the only stops.
Machine time. The full range is 6.87 million candidates: 38 hours on the chain sieve (which in practice would mean days of re-pinning and watching), 23 minutes on the gear sieve, 14 seconds on the ATtiny array. This is the scale of problem the gear machine was built for.
9.6 Rung Five: Diophantine Equations
A Diophantine equation asks for whole-number solutions. Many of the quadratic ones are sieve problems as they stand, because “is this expression a perfect square?” is exactly the test a sieve can pre-screen. kbsriram’s very first firmware for his board (code/01, “The sound of sieving”) solves one with three hard-coded moduli and a buzzer:
y² − 2x² = 4,001.
As a sieve. The candidate is x; the test is whether 4,001 + 2x² is a perfect square. A little reasoning first: 4,001 is odd, so y is odd, so y² ≡ 1 (mod 8); since 4,001 ≡ 1 (mod 8) too, 2x² ≡ 0 (mod 8), and x must be even. The mod-16 ring captures this automatically — it allows exactly the 8 even residues. Running x from 0 to 2,000 (2,001 candidates) and adding moduli in turn, in Python:
Table 5 — As a sieve. The candidate is x; the test is whether 4,001 + 2x² is a perfect square. A little reasoning first: 4,001 is odd, so y is odd, so y² ≡ 1 (mod 8); since 4,001 ≡ 1 (mod 8) too, 2x² ≡ 0 (mod 8), and x must be even. The mod-16 ring captures this automatically — it allows exactly the 8 even residues. Running x from 0 to 2,000 (2,001 candidates) and adding moduli in turn, in Python
| Modulus | Allowed residues | Survivors in 0–2,000 |
|---|---|---|
| 16 | 8 of 16 | 1,001 |
| 9 | 6 of 9 | 667 |
| 5 | 3 of 5 | 400 |
| 7 | 3 of 7 | 171 |
| 11 | 6 of 11 | 95 |
| 13 | 7 of 13 | 52 |
| 17 | 8 of 17 | 22 |
| 19 | 11 of 19 | 14 |
| 23 | 12 of 23 | 9 |
The nine survivors after the 23-ring are x = 38, 52, 280, 350, 500, 602, 1,582, 1,642 and 1,670. Five of them are genuine:
Table 6 — The nine survivors after the 23-ring are x = 38, 52, 280, 350, 500, 602, 1,582, 1,642 and 1,670. Five of them are genuine
| x | y | Check: y² − 2x² |
|---|---|---|
| 38 | 83 | 6,889 − 2,888 = 4,001 |
| 52 | 97 | 9,409 − 5,408 = 4,001 |
| 280 | 401 | 160,801 − 156,800 = 4,001 |
| 350 | 499 | 249,001 − 245,000 = 4,001 |
| 1,642 | 2,323 | 5,396,329 − 5,392,328 = 4,001 |
The other four are false stops that a few more rings would remove. The solutions come in families: multiplying y + x√2 by 3 + 2√2 (the smallest solution of y² − 2x² = 1, the Pell equation) turns one solution into the next. (83 + 38√2)(3 + 2√2) = 401 + 280√2, and (401 + 280√2)(3 + 2√2) = 2,323 + 1,642√2. A sieve does not know about the family structure, but it finds every member in range without being told.
Why the survivors shrink more slowly here. In Rung Three about half the residues were allowed per modulus; here the same is true, but the range is not tuned to the answer, and the equation has infinitely many solutions. The sieve cannot finish with one survivor because there are five genuine answers below 2,000 alone. This is the first hint of the theme of the last section: survivors are made of two parts, genuine solutions (which the sieve must keep) and chance survivors (which more moduli remove).
Machine time. 2,001 candidates: 40 seconds on the chain sieve, under half a second on the gears, 4 milliseconds on the ATtiny array.
9.7 Rung Six: Mersenne and Fermat Numbers
Numbers of the forms 2ⁿ − 1 (Mersenne numbers), 2ⁿ + 1, and 2^(2ⁿ) + 1 (Fermat numbers, Fₙ) were the standard test of every factoring method from Fermat’s day until the computer age. They are attractive for a sieve because their factors are constrained in advance: every prime factor of 2ⁿ − 1, for prime n, has the form 2kn + 1, and is ±1 modulo 8. Those constraints pin down the residues of a in Fermat’s method before the machine ever starts, which shrinks the search enormously.
9.7.1 2⁷⁹ − 1: the gear machine’s showpiece
This is the factorisation that kbsriram’s board reproduces at the end of his project video, and it is laid out in a few lines of Lehmer’s 1933 note. Lehmer began by dividing out the known factor 2,687, leaving
N = (2⁷⁹ − 1) / 2,687 = 224,958,284,260,258,499,201,
a 21-digit number. Its factors must be of the forms 79x + 1 and 8x ± 1. If N = a² − b², those conditions force a into the residue classes a ≡ 159 (mod 6,241 = 79²), a ≡ 0 (mod 3), and a ≡ ±63, ±97, ±159 or ±191 (mod 512). Combining them by the Chinese remainder theorem — Sun Zi’s problem again, on a larger scale — gives a step of 6,241 × 3 × 512 = 9,586,176 and exactly eight starting residues:
a = 9,586,176k + 159, 1,610,337, 4,793,247, 6,403,425, 6,590,655, 7,002,561, 8,987,199 or 9,399,105.
Those eight residues were recomputed for this volume in Python from the three conditions, and they match Lehmer’s list exactly. Lehmer bounded a by 374,930,473,917,097, which gives k ≤ 39,111,579 in each class: “the problem of representing N as the difference of squares was split into 8 parts.” The first two runs found nothing. On the third, “the machine stopped almost at once” at x = 58,088, giving
a = 556,846,584,735, b = 556,644,555,032,
so a − b = 202,029,703 and a + b = 1,113,491,139,767, and
2⁷⁹ − 1 = 2,687 × 202,029,703 × 1,113,491,139,767,
the “13th composite Mersenne number to be completely factored.” Every number in this paragraph was checked in Python, including the primality of all three factors.
Moduli and survivors. All thirty gears of the 1932 machine were used. Following kbsriram’s reconstruction from the 1933 paper, the moduli here are 64, 27, 25, 49 and the 26 primes from 11 to 113. Lehmer’s own 1934 description (gears of 67 to 128 teeth) and the labels on the surviving machine suggest the real set used higher powers such as 128, 81 and 121 (Volume 7), which would only make the sieve stricter. Building each gear’s allowed set from the actual polynomial, as kbsriram’s prob_m79.py does, the 64-, 27- and 79-gears allow every position (their information is already inside the step size), and the other 27 together pass a fraction of about 4.2 × 10⁻⁹ of the candidates. Over a full run of 39.1 million candidates that predicts 0.16 chance stops — so on the run that found the answer, the machine was expected to stop once, at the answer, and it did. Stepping the third class from the start of kbsriram’s range, a Python check of all thirty gears over the first 60,000 candidates finds exactly one survivor.
Machine time. Lehmer’s third run reached the answer after 58,088 candidates: 12 seconds at 5,000 per second, the figure Lehmer published, and 11.3 seconds on kbsriram’s board stepped at the same rate. Counting the two empty runs before it (78.2 million candidates) the real campaign was about 4.3 hours of machine time — consistent with Lehmer’s “the machine ran only a few hours altogether” for the whole month. Covering all eight classes would have taken 17.4 hours on the gears, 72 days on the chain sieve, and 10 minutes on the ATtiny array.
9.7.2 Three more from the same campaign
The 1928 and 1933 papers contain several more problems of the same family, all of which kbsriram reran; the arithmetic below was rechecked in Python.
Table 7 — Three more from the same campaign
| Number | Set-up | Stop | Factors |
|---|---|---|---|
| 2⁹³ + 1, cofactor N = 1,537,228,672,093,301,419 | factors are 93x + 1 and 8x + 1, so a = 69,192k + 1,239,854,518 | k = 6,886 | N = 529,510,939 × 2,903,110,321, so 2⁹³ + 1 = 3² × 529,510,939 × 715,827,883 × 2,903,110,321 |
| 3¹¹¹ + 1, cofactor 2,019,210,335,106,439 | factors are 111n + 1, so a = 49,284x + 44,978,200 | x = 145 | 25,709,599 × 78,539,161 |
| 9,999,000,099,990,001 = (10²⁰ + 1)/(10⁴ + 1) | a = 7,200x + 100,000,601 | x = 400,453 | 1,676,321 × 5,964,848,081 |
The last of these is the one the 1928 paper reports as taking “about 2 hours” on the chain machine; 400,453 candidates at 50 per second is 2 hours 13 minutes, and kbsriram’s rerun took 2 hours 14 minutes. (kbsriram’s driver describes this number as (10¹⁶ + 1)/(10⁴ + 1); that quotient is actually 999,900,009,999. The number in the paper is (10²⁰ + 1)/(10⁴ + 1), which Python confirms divides exactly.) His rerun also stopped once before the answer, at x = 257,749, where a² − N passed all nineteen chains without being a square — a chance survivor the 1928 paper does not mention.
9.7.3 Fermat numbers, and where this rung ends
The Fermat numbers mark the limit of this approach. F₇ = 2¹²⁸ + 1, a 39-digit number, had “withstood many factorization attempts” until Michael Morrison and John Brillhart factored it on an IBM 360/91 at UCLA on the morning of 13 September 1970 (Morrison and Brillhart, “A method of factoring and the factorization of F₇,” Mathematics of Computation 29 (1975) 183–205):
F₇ = 59,649,589,127,497,217 × 5,704,689,200,685,129,054,721,
which Python confirms by multiplication, with both factors testing prime. The smaller factor has 17 digits. In Fermat’s method the two factors are so far apart that a would have to climb from about 1.8 × 10¹⁹ to about 2.9 × 10²¹: some 10²¹ candidates even before the factor constraints are applied, which no sieve of any speed could step. Morrison and Brillhart used a different idea — the continued fraction method, introduced by Lehmer and R. E. Powers in 1931 — which does not search for the answer directly at all, but collects many small partial results and combines them. Their paper is dedicated to D. H. Lehmer on his seventieth birthday. That idea, and the sieve that eventually replaced it, are the subject of Volume 10.
9.8 Rung Seven: Pseudosquares
The last kind of problem is the one the Lehmers’ electronic sieves of the 1960s and 1970s were built for, and it is different in kind from all the others: the sieve is not pre-screening a test, the sieve is the question.
A pseudosquare Lₚ is the smallest positive number that is not a perfect square, is 1 more than a multiple of 8, and is a quadratic residue of every odd prime up to p. It “looks like” a square to every prime up to p, and yet is not one. Pseudosquares matter because they bound how many primes must be checked in certain fast primality tests, and because they measure how long the quadratic-residue pattern of a number can imitate that of a square. D. H. Lehmer, Emma Lehmer and Daniel Shanks tabulated them with the Delay Line Sieve DLS-127 (Lehmer, Lehmer and Shanks, “Integer sequences having prescribed quadratic character,” Mathematics of Computation 24 (1970) 433–451).
As a sieve. Since the condition at 3 forces L ≡ 1 (mod 3), and L ≡ 1 (mod 8), every candidate has the form L = 24K + 1, and the candidate is K. For each prime q from 5 up to p, the ring allows the (q − 1)/2 positions where 24K + 1 is a non-zero square modulo q. No test is run afterwards; a survivor that is not a perfect square is a pseudosquare. The first survivor that is not a square is the answer.
Computing the sequence by brute force in Python, and comparing with OEIS A002189:
Table 8 — Computing the sequence by brute force in Python, and comparing with [OEIS A002189](https://oeis.org/A002189)
| p | Lₚ | p | Lₚ |
|---|---|---|---|
| 3 | 73 | 31 | 515,761 |
| 5 | 241 | 37 | 1,083,289 |
| 7 | 1,009 | 41 | 3,206,641 |
| 11 | 2,641 | 43 | 3,818,929 |
| 13 | 8,089 | 47 | 9,257,329 |
| 17 | 18,001 | 53 | 22,000,801 |
| 19 | 53,881 | 59 | 48,473,881 |
| 23 | 87,481 | 61 | 48,473,881 |
| 29 | 117,049 |
The values agree with the OEIS entry. Note that L₅₉ = L₆₁: the number 48,473,881 happens to be a residue of 61 as well.
Survivors. Here is the instructive part. Sieving K from 0 to 2,019,745 (L up to 48,473,881) with the sixteen primes 5 to 61, the Python run finds 900 survivors. Exactly one is the pseudosquare. The other 899 are perfect squares: a square is a quadratic residue of every prime it is not divisible by, so every square of the form 24K + 1 that avoids the primes 5 to 61 passes every ring. The sieve cannot tell them apart from the answer, and a square check must be run on each one. These are structural survivors — guaranteed by the mathematics, not by bad luck — and no number of additional moduli will ever remove them.
Machine time. 2.02 million candidates: 11 hours on the chain sieve, 7 minutes on the gear sieve, 4 seconds on the ATtiny array. The 1970 problem L₁₂₇ = 196,265,095,009 is another matter: 8.18 billion candidates, which the DLS-127 ran at 1 million per second in about 2.3 hours. kbsriram’s board at 500,000 per second took 4 hours 33 minutes, his log reports; the gear sieve would have needed 19 days and the chain sieve more than five years.
9.9 When Sieving Fails
Every rung so far has been chosen to work. Most problems do not, and the reasons fall into four groups.
9.9.1 Too many candidates
A sieve is a search, and its running time is proportional to the number of candidates. Fermat’s method is efficient only when the two factors are close: the search starts at √N and must reach (p + q)/2. For N = 10,007 × 999,983 (an 11-digit number), that is 404,961 steps — fine. For a 30-digit number whose factors have 11 and 21 digits, (p + q)/2 − √N is about 5 × 10¹⁹ steps, which at the ATtiny array’s 500,000 per second is about 3.2 million years. Lehmer’s tricks — factor constraints, splitting into residue classes, choosing the best quadratic form — can cut the range by large factors, but only by constant factors. The work still grows roughly in proportion to the square root of the number being factored, and that is a losing race against numbers that grow exponentially in their number of digits. This, more than any engineering limit, is why the electromechanical sieves were overtaken. The methods of Volume 10 do not search for the answer; they assemble it from millions of easier pieces.
9.9.2 Too few moduli: chance survivors
Even when the range is affordable, a sieve is only useful if it stops rarely. Each odd prime modulus lets through about half the residues ((p + 1)/2 out of p when the target is a square modulo p). So k moduli pass about one candidate in 2ᵏ, a little more. Every chance survivor costs a human — or a host computer — a full check.
The chart makes a point about the history of the machines. A faster sieve needs more moduli, not fewer, because it meets chance survivors faster. Computed from the actual allowed links of two 1928 problems (the 9,999,000,099,990,001 run and the second sum-of-two-squares run), the nineteen chains passed about one candidate in 220,000 to 260,000 — at 50 per second, a chance stop somewhat less than once an hour, which an operator watching the machine could live with. Thirty gears at 5,000 per second, on the bare-prime estimate in the chart, give about one chance stop every eight hours. kbsriram’s array has thirty cards and so carries the gear machine’s filtering power at a hundred times its speed: about twelve chance stops an hour on a problem with no extra structure. That is nothing for a host computer that checks a survivor in microseconds, and it would have swamped a human.
Two of kbsriram’s reruns show this in miniature. The 1928 sum-of-two-squares run and the 9,999,000,099,990,001 run both produced a survivor that the nineteen 1928 chains could not reject: in the first case N − a² was a residue of every prime up to 67 and failed only at 73. The cure is always the same — another modulus — and the 1928 machine had run out of them.
9.9.3 Moduli that do nothing
Some moduli carry no information for a particular problem. In the 2⁷⁹ − 1 run the 64-, 27- and 79-gears allowed every position, because the step 9,586,176 = 2⁹ × 3 × 79² had already used up all that those moduli could say. In the 1928 runs the 64- and 25-chains were pinned in every link for the same reason. A machine with thirty moduli may have only twenty-seven working on a given problem, and its filtering power is computed from those alone.
9.9.4 Survivors that are not chance
The pseudosquare search showed the worst case: survivors that pass every modulus because the question cannot distinguish them from the answer. Every perfect square passes a quadratic-residue sieve. Trivial solutions of an equation (x = 0, y = √N) pass a sieve for that equation. In factoring by quadratic forms, the representation that comes from N itself, or from a known factor, will always appear. These have to be predicted and removed by reasoning before the run, or recognised and skipped during it. kbsriram’s drivers include a certain_survivors function for exactly this reason; for the pseudosquare problems it is the count of squares that must be expected.
9.9.5 Summary of the ladder
Table 9 — Summary of the ladder
| Problem | Candidate | Moduli used | Survivors | Chain (50/s) | Gears (5,000/s) | ATtiny (500k/s) |
|---|---|---|---|---|---|---|
| Primes to 100 | n, 2–100 | 4 (exclusive) | 25 primes | 2 s | 20 ms | 0.2 ms |
| Sun Zi 3/5/7 | n, 0–104 | 3 | 1 per 105 | 2 s | 21 ms | 0.2 ms |
| Fermat, 10,012,997 | a, 337 values | 7 | 1 | 7 s | 67 ms | 0.7 ms |
| y² − 2x² = 4,001 | x, 0–2,000 | 9 | 9 (5 genuine) | 40 s | 0.4 s | 4 ms |
| 3¹¹¹ + 1 cofactor | x, to 145 | 19 | 1 | 3 s | 29 ms | 0.3 ms |
| 2⁹³ + 1 cofactor | k, to 6,886 | 30 | 1 | 2 min | 1.4 s | 14 ms |
| 2⁷⁹ − 1, third run | x, to 58,088 | 30 (27 active) | 1 | 19 min | 12 s | 0.1 s |
| 4,500,383,424,373 = a² + b² | x, to 157,044 | 19 | 2 (1 genuine) | 52 min | 31 s | 0.3 s |
| 9,999,000,099,990,001 | x, to 400,453 | 19 | 2 (1 genuine) | 2.2 h | 80 s | 0.8 s |
| L₅₉ = 48,473,881 | K, 2.02 million | 16 | 900 (1 non-square) | 11 h | 7 min | 4 s |
| x² + 7y² for M, full range | k, 6.87 million | 30 | 2 | 38 h | 23 min | 14 s |
| 2⁷⁹ − 1, all eight classes | 313 million | 30 | — | 72 days | 17 h | 10 min |
| L₁₂₇ = 196,265,095,009 | K, 8.18 billion | 29 | — | 5.2 years | 19 days | 4.5 h |
The table is the argument of the whole series in one place. Each machine was a hundred times faster than the last, and each opened a band of problems the previous one could not touch. None of them changed the shape of the curve: double the digits of N in a factoring problem and a direct search becomes impossibly long on any of them. The next volume follows the line of thought that escaped that curve — and finds, at the heart of the fastest factoring methods ever devised, the same act of throwing numbers away.
Next — Volume 10: Sieves in Cryptography — From Lehmer to the Number Field Sieve. Why factoring became a matter of security in 1977, and how the same idea of throwing numbers away grew into the quadratic sieve, the number field sieve, Shamir’s optical TWINKLE, the RSA-260 record of 2026, and the post-quantum standards that end the story.
Comments (0)