# [Haskell-cafe] Even better Eratosthenes sieve and lucky numbers

oleg at pobox.com oleg at pobox.com
Mon Feb 12 21:37:46 EST 2007

We further simplify the previously posted genuine sieve algorithm and
generalize it to the finding of lucky numbers.

We observe that we only need to store marks _signifying_ the integers,
but never the integers themselves. Thus we arrive at the algorithm
that is distinguished from all previously posted by:
(i) doing no operations on composite numbers
(ii) using neither multiplication nor division nor the
remainder operation
(iii) using neither general addition nor the comparison.

The algorithm only relies on the successor, predecessor and zero
comparison. The predecessor can be easily eliminated. Thus the
algorithm can be used with Church and Peano numerals, or members of
Elliptic rings, where zero comparison and successor take constant
time but other arithmetic operations are more involved.

> -- repl_every_n n l replaces every (n+1)-th element in a list (_:l)
> -- with False
> repl_every_n :: Int -> [Bool] -> [Bool]
> repl_every_n n l = repl_every_n' n l
>  where repl_every_n' 0 (_:t) = False: repl_every_n n t
>        repl_every_n' i (h:t) = h:     repl_every_n' (pred i) t
>
> primes = 2:(loop 3 (repeat True))
>  where loop n (False:t) = loop (succ (succ n)) t
>        loop n (_:t)  = n:(loop (succ (succ n)) (repl_every_n (pred n) t))
>
> main = putStrLn \$ "Last 10 primes less than 10000: " ++
>        show (take 10 \$ reverse  \$ takeWhile (< 10000) primes)

Last 10 primes less than 10000:
[9973,9967,9949,9941,9931,9929,9923,9907,9901,9887]

The algorithm easily generalizes to lucky numbers
http://mathworld.wolfram.com/LuckyNumber.html
http://www.research.att.com/~njas/sequences/A000959

> -- Consider the series of odd numbers starting with k=1; the k-th number
> -- is 2k-1.
> -- If we start the elimination from the very beginning of the series (k=1),
> -- we start with the phase (2k-2) and eliminate with the step 2k-2.
> -- If we start at the number k+1, the phase becomes (2k-2-k) = (k-2).
> -- Thus the starting phase gets incremented by 1 as we move up the series.
> -- However, already eliminated numbers don't count, so as we skip
> -- over eliminated numbers, we increment the phase by 2
> lucky = 1:(loop 3 0 (repeat True))
>  where loop n i (False:t) = loop (succ (succ n)) (succ (succ i)) t
>        loop n i (_:t)  = n:(loop (succ (succ n)) (succ i)
>                             (repl_every_np i (pred n) t))
>
> -- repl_every_np i n eliminates (marks as False) every (n+1)-th number
> -- starting with the phase i.
> -- Already eliminated numbers don't count.
> repl_every_np :: Int -> Int -> [Bool] -> [Bool]
> repl_every_np i n (False:t) = False: repl_every_np i n t
> repl_every_np 0 n (_:t) = False: repl_every_np n n t
> repl_every_np i n (_:t) = True:  repl_every_np (pred i) n t

*Main> take 10 lucky
[1,3,7,9,13,15,21,25,31,33]
*Main> takeWhile (<100) lucky
[1,3,7,9,13,15,21,25,31,33,37,43,49,51,63,67,69,73,75,79,87,93,99]