99 questions/Solutions/31

From HaskellWiki
Jump to navigation Jump to search
The printable version is no longer supported and may have rendering errors. Please update your browser bookmarks and please use the default browser print function instead.

(**) Determine whether a given integer number is prime.

isPrime :: Integral a => a -> Bool
isPrime p = p > 1 && 
            (all ((/= 0).(p `rem`)) $ candidateFactors p)

candidateFactors p = takeWhile ((<= p).(^2)) [2..]

Well, a natural number p is a prime number if it is larger than 1 and no natural number n >= 2 with n^2 <= p is a divisor of p. That's exactly what is implemented: we take the list of all integral numbers starting with 2 as long as their square is at most p and check that for all these n there is a non-zero remainder concerning the division of p by n.

However, we don't actually need to check all natural numbers <= sqrt P. We need only check the primes <= sqrt P:

{-# OPTIONS_GHC -O2 -fno-cse #-}
candidateFactors p = let z = floor $ sqrt $ fromIntegral p + 1 in
                     takeWhile (<= z) primesTME

-- tree-merging Eratosthenes sieve
--  producing infinite list of all prime numbers
primesTME = 2 : gaps 3 (join [[p*p,p*p+2*p..] | p <- primes'])
  where
    primes' = 3 : gaps 5 (join [[p*p,p*p+2*p..] | p <- primes'])
    join  ((x:xs):t)        = x : union xs (join (pairs t))
    pairs ((x:xs):ys:t)     = (x : union xs ys) : pairs t
    gaps k xs@(x:t) | k==x  = gaps (k+2) t 
                    | True  = k : gaps (k+2) xs

-- duplicates-removing union of two ordered increasing lists
union (x:xs) (y:ys) = case (compare x y) of 
           LT -> x : union  xs  (y:ys)
           EQ -> x : union  xs     ys 
           GT -> y : union (x:xs)  ys

The tree-merging Eratosthenes sieve here seems to strike a good balance between efficiency and brevity. More at Prime numbers haskellwiki page.