Factor Tables
May 8, 2012
Before the dawn of computers, most computations were done with the aid of tables: logarithm tables, sine tables, and so on. These tables were ubiquitous, indispensable, and riddled with errors. Number theorists who needed to factor numbers used tables of the least prime factor of a number. The oldest such table dates to 1603 (it contained the least prime factor of all numbers to 750), and new tables were being constructed as late as Derrick N. Lehmer’s table of least prime factors to ten million in 1909 (he was the father of Derrick H. Lehmer); Maarten Bullynck gives the history. Here’s a sample page from a large table, showing the least prime factors of all numbers less than a thousand; numbers divisible by 2 and 5 are omitted, and primes are skipped:
0 1 2 3 4 5 6 7 8 9
1 -- 3 7 -- 3 -- -- 3 17
3 -- -- 7 3 13 -- 3 19 11 3
7 -- -- 3 -- 11 3 -- 7 3 --
9 3 -- 11 3 -- -- 3 -- -- 3
11 -- 3 -- -- 3 7 13 3 -- --
13 -- -- 3 -- 7 3 -- 23 3 11
17 -- 3 7 -- 3 11 -- 3 19 7
19 -- 7 3 11 -- 3 -- -- 3 --
21 3 11 13 3 -- -- 3 7 -- 3
23 -- 3 -- 17 3 -- 7 3 -- 13
27 3 -- -- 3 7 17 3 -- -- 3
29 -- 3 -- 7 3 23 17 3 -- --
31 -- -- 3 -- -- 3 -- 17 3 7
33 3 7 -- 3 -- 13 3 -- 7 3
37 -- -- 3 -- 19 3 7 11 3 --
39 3 -- -- 3 -- 7 3 -- -- 3
41 -- 3 -- 11 3 -- -- 3 29 --
43 -- 11 3 7 -- 3 -- -- 3 23
47 -- 3 13 -- 3 -- -- 3 7 --
49 7 -- 3 -- -- 3 11 7 3 13
51 3 -- -- 3 11 19 3 -- 23 3
53 -- 3 11 -- 3 7 -- 3 -- --
57 3 -- -- 3 -- -- 3 -- -- 3
59 -- 3 7 -- 3 13 -- 3 -- 7
61 -- 7 3 19 -- 3 -- -- 3 31
63 3 -- -- 3 -- -- 3 7 -- 3
67 -- -- 3 -- -- 3 23 13 3 --
69 3 13 -- 3 7 -- 3 -- 11 3
71 -- 3 -- 7 3 -- 11 3 13 --
73 -- -- 3 -- 11 3 -- -- 3 7
77 7 3 -- 13 3 -- -- 3 -- --
79 -- -- 3 -- -- 3 7 19 3 11
81 3 -- -- 3 13 7 3 11 -- 3
83 -- 3 -- -- 3 11 -- 3 -- --
87 3 11 7 3 -- -- 3 -- -- 3
89 -- 3 17 -- 3 19 13 3 7 23
91 7 -- 3 17 -- 3 -- 7 3 --
93 3 -- -- 3 17 -- 3 13 19 3
97 -- -- 3 -- 7 3 17 -- 3 --
99 3 -- 13 3 -- -- 3 17 29 3
For example the table shows that the least prime factor of 923, in the column headed 9 and row headed 23, is 13, and 997 is the greatest prime less than a thousand. To factor a number, find its least prime factor in the table, compute the remaining cofactor by division, and repeat until the cofactor is prime.
When building the table, the least prime factor is computed by sieving, not by trial division. The setup is the same as the Sieve of Eratosthenes, except that integers are used instead of booleans, and each item in the sieve is initialized to 1. Then each successively-smallest prime p is sieved, but instead of changing true to false, each 1 in the chain of multiples of p is changed to p (other values are ignored). The same optimizations as the normal sieve — odd numbers only, start at the square of p, and stop when p2 is greater than n — apply here.
Your task is to write a function that sieves for least prime factors, and to use that function to write a program that builds factor tables as illustrated above. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.
Programming Praxis 
Forth implementation:
1000000 constant LPF_SIZE : sq ( x -- x^2 ) s" dup *" evaluate ; immediate : create-table ( compile: sz -- ; run: n -- addr ) create here over cells allot swap 0 DO 1 !+ LOOP drop does> swap 1- cells + ; LPF_SIZE create-table lpf : sieve-lpf ( -- ) 2 BEGIN dup lpf @ 1 = IF dup sq BEGIN dup LPF_SIZE <= WHILE dup lpf @ 1 = IF 2dup lpf ! THEN over + REPEAT drop THEN 1+ dup sq LPF_SIZE > UNTIL drop ; : .factors ( n -- ) BEGIN dup lpf @ dup 1 <> WHILE dup . / REPEAT drop . ; : .header ( n00 k -- ) cr ." " 0 DO dup 5 .r 1+ LOOP drop cr ; : show? ( n -- ? ) dup 1 and 0<> swap 5 mod 0<> and ; : .row ( n00 k row -- ) dup 2 .r -rot bounds DO i 100 * over + lpf @ dup 1 = IF ." --" drop ELSE 5 .r THEN LOOP drop cr ; : .table ( n00 k -- ) 2dup .header 100 0 DO i show? IF 2dup i .row THEN LOOP 2drop ; sieve-lpfTesting:
Python implementation:
Testing:
Python implementation: http://pastebin.com/CFP3S2HV
Testing: http://pastebin.com/T1Q0W8pP
c#: http://pastebin.com/41KxNmpP
I didn’t notice that numbers divisible by 2 and 5 are omitted, and I replaced — with blanks
Another Python implementation.
pagenr = 1000 # 1000 numbers on a page; 1: 0-1000 N = 1000 last = pagenr * N lowprime = ["--"] * (last + 1) for i in range(3, last + 1, 2): if lowprime[i] == "--": for j in xrange(min(i * i, last), last + 1, 2 * i): if lowprime[j] == "--": lowprime[j] = i fs = lambda x: "{0:>6s}".format(str(x)) print fs(""), for k in xrange((last - N) // 100, (last) // 100): print fs(k), for row in xrange(1, 100, 2): if row % 5: print print fs(row), for k in lowprime[row + last - N::100]: print fs(k),