99 questions/Solutions/39

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.

(*) A list of prime numbers.

Given a range of integers by its lower and upper limit, construct a list of all prime numbers in that range.

Solution 1:

primesR :: Integral a => a -> a -> [a]
primesR a b = filter isPrime [a..b]

If we are challenged to give all primes in the range between a and b we simply take all number from a up to b and filter the primes out.

Solution 2:

primes :: Integral a => [a]
primes = let sieve (n:ns) = n:sieve [ m | m <- ns, m `mod` n /= 0 ] 
         in sieve [2..]

primesR :: Integral a => a -> a -> [a]
primesR a b = takeWhile (<= b) $ dropWhile (< a) primes

Another way to compute the claimed list is done by using the Sieve of Eratosthenes. The primes function generates a list of all (!) prime numbers using this algorithm and primesR filter the relevant range out. [But this way is very slow and I only presented it because I wanted to show how nicely the Sieve of Eratosthenes can be implemented in Haskell :)]

Solution 3:

Use the proper Sieve of Eratosthenes from e.g. 31st question's solution (instead of the above sieve of Turner), adjusted to start its multiples production from the given start point:

{-# OPTIONS_GHC -O2 -fno-cse #-}
-- tree-merging Eratosthenes sieve, primesTME of haskellwiki/prime_numbers, 
--  adjusted to produce primes in a given range
primesR a b 
  | b<a || b<2 = []
  | otherwise  = 
     (if a <= 2 then [2] else []) ++
     gaps a' (join [[x,x+step..b] | p <- takeWhile (<= z) (tail primesTME)
                    , let q = p*p ; step = 2*p
                          x = if a' <= q then q else snapUp a' q step ])
  where
    a'      = if a<=3 then 3 else (if even a then a+1 else a)
    z       = floor $ sqrt $ fromIntegral b + 1
    join  (xs:t)    = union xs (join (pairs t))
    join  []        = []
    pairs (xs:ys:t) = (union xs ys) : pairs t
    pairs  t        = t
    gaps k xs@(x:t) | k==x  = gaps (k+2) t 
                    | True  = k : gaps (k+2) xs
    gaps k []       = [k,k+2..b]

    snapUp v origin step = let r = rem (v-origin) step
                           in if r==0 then v else v-r+step
    -- 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
    union  a      b     = a ++ b

(This turned out to be quite a project, with a few quite subtle points.) It should be much better then taking a slice of a full sequential list of primes, as it won't try to generate any primes between the square root of b and a. To wit,

> primesR 10100 10200                                            -- Sol.3
[10103,10111,10133,10139,10141,10151,10159,10163,10169,10177,10181,10193]
(6,038 reductions, 11,923 cells)

> takeWhile (<= 10200) $ dropWhile (< 10100) $ primesTME  -- TME: of Q.31
[10103,10111,10133,10139,10141,10151,10159,10163,10169,10177,10181,10193]
(140,313 reductions, 381,058 cells)

> takeWhile (<= 10200) $ dropWhile (< 10100) $ sieve [2..]       -- Sol.2
     where sieve (n:ns) = n:sieve [ m | m <- ns, m `mod` n /= 0 ]
[10103,10111,10133,10139,10141,10151,10159,10163,10169,10177,10181,10193]
(54,893,566 reductions, 79,935,263 cells, 6 garbage collections)

> filter isPrime [10100..10200]                                  -- Sol.1
[10103,10111,10133,10139,10141,10151,10159,10163,10169,10177,10181,10193]
(34,694 reductions, 55,146 cells)  -- isPrime: Q.31-Sol.1 using primesTME
(15,750 reductions, 29,292 cells)  -- isPrime: Q.31-Sol.2 using primesTME

(testing with Hugs of Nov 2002).