[Haskell-cafe] A tale of Project Euler

[Sebastian Sylvan]

On Nov 27, 2007 8:54 PM, Olivier Boudry <olivier.boudry at gmail.com> wrote:
>

> Hi Andrew,
>
>  I don't remember who I stole this prime computation from, but it is very
> fast (10001's prime in 0.06 sec with Int and 0.2 with Integer on my machine)
> and not overly complex:
>
>  primes :: [Integer]
>  primes = 2 : filter (l1 . primeFactors) [3,5..]
>      where
>          l1 (_:[]) = True
>          l1      _ = False
>
>  primeFactors :: Integer -> [Integer]
>  primeFactors n = factor n primes
>      where
>          factor _ [] = []
>          factor m (p:ps) | p*p > m        = [m]
>                          | m `mod` p == 0 = p : factor (m `div` p) (p:ps)
>                          | otherwise      = factor m ps
>
> I used it a lot while playing with Euler Project. Many of the problems
> require prime calculation.
>


That is indeed a nice and clear version that's pretty fast. It's
basically the same as the C version except "backwards" (i.e. examine a
number and work backwards through its divisors, rather than filling in
a "map" of all multiples of a known prime).
Let me suggest the following slight modification (primeFactors in your
version doesn't actually return prime factors - it returns prime
factors *or* a list of just the number itself),

primes :: [Integer]
primes = 2 : filter (null . primeFactors) [3,5..]

primeFactors :: Integer-> [Integer]
primeFactors n = factor n primes
    where
        factor m (p:ps) | p*p > m        = []
                        | m `mod` p == 0 = p : factor (m `div` p) (p:ps)
                        | otherwise      = factor m ps


This is roughly 35% faster on my machine with GHC 6.7.20070730 too,
but the point wasn't to make it faster, but clearer.
-- 
Sebastian Sylvan
+44(0)7857-300802
UIN: 44640862