Enumerating primes on a Turing machine in nearly linear time, unconditionally

· 12 min · 2401 words · nor

TL;DR

Harvey [1] recently gave a conditional algorithm that lists the primes below NN in N(log⁡log⁡N)1+o(1)N (\log \log N)^{1+o(1)} bit operations on a multitape Turing machine, assuming a conjecture about crowded short intervals (without it, the algorithm may report failure, though it never returns a wrong list). Unconditionally, his deterministic algorithm costs about log⁡N\sqrt{\log N} more. We show an algorithm with N(log⁡log⁡N)1+o(1)N (\log \log N)^{1+o(1)} deterministically and with no conjecture (see how this was done).

There are two documents: an intuitive explanation understandable by someone outside this field, which is the place to start, and the paper, which has the full proof (and an overview of the proof in the first section). This post is about the main idea, and some things around it.

Some context

The sieve of Eratosthenes (3rd-century BCE) is remarkably fast for an algorithm that ancient. With random access to memory and word size logarithmic in input (both conditions are important), it takes O(Nlog⁡log⁡N)O(N \log \log N) operations, and the best known bound stayed at O(N/log⁡log⁡N)O(N / \log \log N) from Pritchard's 1981 paper until this year - Pritchard even conjectured that it was optimal.

On a multitape Turing machine (a much more conservative setting), where each tape head moves at most one cell per operation, so accessing distant memory locations can require many operations, things change a lot.

Naively implementing the sieve of Eratosthenes gives a bound of O(N3/2/log⁡N)O(N^{3/2} / \log N), which is nowhere near close to what we want. The best known algorithm in this setting was Sergeev's O(Nlog⁡N)O(N \log N) bit operations.

Harvey's paper [1] beat both Sergeev's algorithm and (when converted to RAM-machine bounds) Pritchard's algorithm, by a lot: his conditional bound is about log⁡N(log⁡log⁡N)1+o(1)\log N \, (\log \log N)^{1+o(1)} bit operations per prime (and there are around N/ln⁡NN / \ln N of them below NN), barely more than a single multiplication of two log⁡N\log N-bit integers (he also has a randomized Las Vegas version of the algorithm which adds roughly an extra log⁡log⁡N\log \log N, and a deterministic unconditional version that has roughly an extra log⁡N\sqrt{\log N}). Just printing out a list of all the primes below NN takes linear time in NN, so the conditional version (and the randomized version) is quite close to optimality.

It can be shown that, by using the logarithmic word size to help simulate the multitape Turing machine faster on a RAM machine, we can shave off a log⁡N\log N factor, which, like Pritchard's algorithm, makes [1]'s algorithms faster than even reading the first NN natural numbers (and faster than Pritchard's algorithm, by a positive power of log⁡N\log N, without doing anything RAM-machine-specific other than this one reduction).

Unfortunately, the best result in [1] is conditional on a conjecture that seems out of reach of current methods to prove.

Our result removes this dependence on the conjecture, and answers one of the open questions in his paper: whether N(log⁡N)o(1)N(\log N)^{o(1)} is possible deterministically and unconditionally. In fact, we obtain the stronger bound N(log⁡log⁡N)1+o(1)N(\log\log N)^{1+o(1)}.

The idea behind the approach

Harvey's core algorithm is a kind of compressed sensing machine. It finds the odd square-primes m2pm^2 p (which can be used to get primes fast enough, by another result in [1]) block by block, and stores each block in a "sketch" (summary) with some capacity, which has to be about a constant times the average number of square-primes in a block. If a block fits, the sketch decodes correctly to the odd square-primes in that range. If it overflows, the sketch can decode to something incorrect, and we would not know when that happens (unless we did an expensive test).

If square-primes behaved like random points, a block would overflow with probability exponentially small in its mean occupancy, so essentially never - and Harvey's conjecture is roughly about this, which seems far out of reach of current methods. But we don't necessarily need "never". [1] also has an expensive direct test, and we can afford to run it on about a (log⁡N)−12(\log N)^{-12} fraction of the numbers. So being able to say something about almost all blocks would be enough - which we can actually hope to prove - as long as the algorithm can also find the bad blocks.

We never count square-primes in a block. Instead, we bound them by a sieve weight - the sieve of Eratosthenes, truncated carefully enough that it stays an upper bound, using what is called the beta sieve - which has the same density as the primes up to a constant, and is a divisor sum over a small list of divisors. Essentially, the idea is to show that there is a nonnegative weight that is always greater than the indicator function of the square-primes we care about (with some small caveats), and then we sum it over blocks/intervals of interest. Then there are two questions: is this weight itself rarely large, and can we compute it quickly? The answer to both is "not exactly, but on average", and it turns out that's basically enough to give us something useful to work with.

Figure 1: Left: sieving by 2, 3 and 5 keeps the primes, along with composites like 49 and 77 whose prime factors are all larger - an upper bound for the primes, which is what the weight is. Right: written as inclusion-exclusion (here for primes below 30), simply dropping the terms with divisors above 1000 can push the weight below the true value, even negative, for about 1% of n, while the beta-sieve (here with β = 2) never goes below it.
Figure 1: Left: sieving by 2, 3 and 5 keeps the primes, along with composites like 49 and 77 whose prime factors are all larger - an upper bound for the primes, which is what the weight is. Right: written as inclusion-exclusion (here for primes below 30), simply dropping the terms with divisors above 1000 can push the weight below the true value, even negative, for about 1% of n, while the beta-sieve (here with β = 2) never goes below it.

Rarely large. The mean of the weight already fits comfortably in the capacity, so only fluctuations can hurt. And there's a free parameter that we can vary: the core doesn't care how long blocks are, as long as the capacity scales along and the block sizes are at most polylogarithmic in NN. Longer blocks average out better - the fluctuations of a block sum grow like the square root of its length (up to certain other logarithmic factors), while overflowing needs a deviation proportional to its length - so a variance bound and Chebyshev are enough. The variance is small because the counting errors for different divisors have limited correlations. For coprime divisors, the Chinese remainder theorem makes these correlations exactly zero over complete periods, and the remaining correlations can also be bounded.

Figure 2: Odd square-primes in blocks of length T across \([2^{27}, 2^{28})\), relative to the mean, against random points with the same mean, and an illustrative capacity of \(1.5 \times\) the mean. Longer blocks concentrate, so overflowing gets rarer. The real counts are even more concentrated than random points (their variance is only about 0.5-0.6 times the mean).
Figure 2: Odd square-primes in blocks of length T across [227,228)[2^{27}, 2^{28}), relative to the mean, against random points with the same mean, and an illustrative capacity of 1.5×1.5 \times the mean. Longer blocks concentrate, so overflowing gets rarer. The real counts are even more concentrated than random points (their variance is only about 0.5-0.6 times the mean).

Computing it. Over a block, it turns out that the weight counts pairs (d,k)(d, k) with dkdk in the block (where dd and kk can vary sufficiently independently). There are far too many pairs to list, but log⁡dk=log⁡d+log⁡k\log dk = \log d + \log k, so this becomes a convolution in log-space, which we can compute using a constant number of big-integer multiplications. Since this uses finite precision, we enclose product locations in intervals and treat positive and negative contributions conservatively, so the resulting block counts are always upper bounds. Rounding errors only happen near block boundaries, which are a vanishing fraction of the line (if we choose the block size carefully), so the total slack is small even though some blocks get terrible bounds. Only the total matters, since overshoots being nonnegative implies that only a few blocks can have a big overshoot (this is like a Markov bound).

Figure 3: Each output of the convolution only knows where its products are up to a short interval. Inside a block, its net count goes to that block. If it crosses a boundary, the positive count goes to both blocks and the negative count to neither, so the bound never undercounts, and all the slack sits near block boundaries.
Figure 3: Each output of the convolution only knows where its products are up to a short interval. Inside a block, its net count goes to that block. If it crosses a boundary, the positive count goes to both blocks and the negative count to neither, so the bound never undercounts, and all the slack sits near block boundaries.

Everything else is book-keeping and choosing powers of log⁡N\log N so that all of this fits, including a small arithmetic trick about why the block length should be a power of 22. The intuitive explanation goes through all of it, and is written in a top-down manner.

Figure 4: The whole algorithm.
Figure 4: The whole algorithm.

Some interesting observations

"Almost all" is enough, if you can find the exceptions. In [1], the repair step helps us go from "almost all" results to worst-case guarantees. The only extra requirement is that the exceptions can be located cheaply, which is where most of our work goes. The bound we end up needing is also not super exotic - mean-square estimates of this kind are classical, e.g., a closely related estimate appears in Matomäki's work [5] on almost primes in almost all very short intervals.

Exact but dense, or approximate but sparse. [1]'s unconditional algorithm also uses a sieve-based proxy, but computes it exactly, which only works with a small sieving limit, so the proxy is too dense, which makes the algorithm slower. Our algorithm sieves up to a power of NN, so the density is the right order, and it only gets computed approximately - which is not an issue, since only the total error matters (and we bound it to be small).

Figure 5: A small-scale version of the majorant (primes below 30, divisors up to 20000, \(\beta\) = 2) against the true count, on blocks of length \(2^{12}\) in \([2^{23}, 2^{24})\). It never drops below the count, and stays within a constant factor of it - between \(2.3 \times\) and \(3 \times\) on every block.
Figure 5: A small-scale version of the majorant (primes below 30, divisors up to 20000, β\beta = 2) against the true count, on blocks of length 2122^{12} in [223,224)[2^{23}, 2^{24}). It never drops below the count, and stays within a constant factor of it - between 2.3×2.3 \times and 3×3 \times on every block.

The bottleneck to improving complexity is now the core algorithm. Everything we did is linear time. [1] gives a heuristic argument that Nlog⁡log⁡NN \log \log N is about the best possible for this kind of approach: the core recovers about Nlog⁡log⁡N/log⁡NN \log \log N / \log N bits of information, and polynomial multiplication should cost a log⁡N\log N factor over that. So for this approach, the remaining gap is only the (log⁡log⁡N)o(1)(\log \log N)^{o(1)}. This log⁡N\log N factor might need to be reconsidered in case methods being used in recent proposed improvements to integer multiplication help.

The constants are very large, as usual. This is a galactic algorithm. Blocks have length about (log⁡N)18(\log N)^{18}, so for N=10100N = 10^{100}, a single block would have about 104510^{45} numbers (even if the constant in front of the power was 1), and the constant in the capacity is 400000400000 times another large constant. Similarly, the probabilistic algorithm implicitly assumes N>exp⁡(80e3)≈10698N > \exp(80 e^3) \approx 10^{698}. None of this is practical, and none of it has been optimized (in particular, we use a constant of 1/10001/1000 for something where something as large as 1/4−ε1/4 - \varepsilon might potentially work, but it was easier to cite the result from Matomäki with the same constant). This is also generally how results in the field look like, so no one's really at fault.

Problems that are still open

  • Removing the (log⁡log⁡N)o(1)(\log \log N)^{o(1)} factor. Harvey's recent note [4] on the proposed o(nlog⁡n)o(n \log n) integer multiplication algorithm speculates that a (conditional) faster multiplication algorithm might help here, after going through the whole complexity analysis.
  • A version that uses much less than linear memory.
  • The estimate for crowded blocks seems to be enough to bring Harvey's Las Vegas algorithm down to the same expected time as well, keeping the rest of his randomized algorithm, but we haven't written this up.

How this was done

[1] did most of the heavy lifting - the core algorithm, the direct test and the conversion to primes are all Harvey's. The main proof was mostly by GPT-6 Astra, heavily steered until I could understand everything, and Claude Opus 5.5 helped turn my notes into the intuitive explanation (which I had to iterate upon) and draft this post. The LLM-written proofs were, as usual, hard to read, so most of the few days this took went into verification and exposition (the latter taking more time than the former). I've verified both documents manually as well as using LLMs, but there is no Lean formalization, so if you find a mistake, please let me know.

References

  1. D. Harvey, Faster enumeration of primes, arXiv:2606.22851, 2026.
  2. Enumerating primes in nearly linear time, intuitively, notes, 24 October 2026.
  3. Deterministic enumeration of primes in N(log⁡log⁡N)1+o(1)N (\log \log N)^{1+o(1)} bit operations, manuscript, 24 October 2026.
  4. D. Harvey, On the recently proposed o(nlog⁡n)o(n \log n) integer multiplication algorithm, note, 20 October 2026.
  5. K. Matomäki, Almost primes in almost all very short intervals, J. London Math. Soc. 106 (2022), 1061-1097.

Cite this post

@online{prime-enumeration,
  author = {nor},
  title = {Enumerating primes on a Turing machine in nearly linear time, unconditionally},
  year = {2026},
  month = {10},
  day = {27},
  url = {https://nor-blog.pages.dev/posts/2026-10-27-prime-enumeration/},
}