All pairs of consecutive 47-smooth numbers, calculated by my Python code for Lehmer's improvement to Størmer's theorem. That is, this file lists pairs x, x+1 such that all prime factors of both x and x+1 are at most 47. There are 1502 such pairs, the largest of which is formed by 1109496723125 = 5⁴ × 7 × 17 × 19² × 31² × 43 and 1109496723126 = 2 × 3 × 11² × 23 × 29² × 41² × 47.

It took nearly two days on a dual-core Pentium, so the next prime 53 may be in reach but would likely take more computer time than I'm willing to devote to it. Judging by A002072, faster codes are out there.

Størmer's theorem is really not so much a theorem but an algorithm, for finding pairs of consecutive smooth numbers. I implemented Lehmer's improved algorithm, in Python, as I needed an implementation to generate some examples while expanding the WP article (and maybe more importantly, to check my understanding of Lehmer's paper). It includes some careful bignum implementation of rounding quadratic irrationals, finding continued fraction expansions of them, and solving Pell's equation, which might be useful for other purposes.

I just merged today two articles in Wikipedia on the same problem, regular numbers, that is, numbers of the form 2ⁱ3^j5^k. The ancient Babylonians liked these numbers because their reciprocals had terminating sexagesimal representations; more recently, "Hamming's problem" is the task of generating these numbers efficiently. Dijkstra described a very simple but inefficient functional solution: if H is the sequence of numbers to be generated, then H consists of the number 1 concatenated with the merger of 2H, 3H, and 5H. But there are two flaws with this. First, it generates the same numbers many times: e.g. 6 is generated both as 2*3 and as 3*2. This can be alleviated by letting A be the powers of 2, B be A merged with 3B, and H be B merged with 5H. But even if that improvement is made, and one does the merges in the most obvious and natural lazy-functional way, the sequence is generated inefficiently: the number of merges a generated number is involved in is O(j+k), and all these merges cause generating n numbers in this way to take Θ(n^4/3) arithmetic operations.

Below is an approach that instead generates n numbers in a linear number of operations. The basic principle is, when solving equations like H = merge(B, 5*H), instead of generating 5*H using a recursive instance of the same algorithm, to cache a single list of the whole sequence of H and reuse that stored sequence when computing 5*H.

( Cut for Python code )

However, the pure functional approach uses less space to generate the same sequence, O(n^2/3) instead of O(n) for the caching approach. I can find algorithms that generate the sequence in O(n^2/3) space and O(n log n) time (essentially, by turning it into a graph shortest path problem in an infinite grid graph, applying Dijkstra's algorithm for shortest paths, and only keeping track of nodes on the search frontier), but is it possible to solve optimally both in time and space?

ETA: Yes, it is, by a very simple change to the above code. Simply throw away the parts of the cached sequence that the back pointer has moved past. ( Read more... )

ETA2, March 14: Writing in a more imperative style actually gives tighter code, though perhaps it is not as readable: ( Read more... )