# 99 questions/Solutions/31

### From HaskellWiki

< 99 questions | Solutions(Difference between revisions)

(putting in a better alternative for allPrimes) |
|||

Line 4: | Line 4: | ||

isPrime :: Integral a => a -> Bool |
isPrime :: Integral a => a -> Bool |
||

isPrime p = p > 1 && |
isPrime p = p > 1 && |
||

− | (all (\n -> p `mod` n /= 0 ) |
+ | (all ((/= 0).(p `rem`)) $ candidateFactors p) |

− | $ takeWhile (\n -> n*n <= p) [2..]) |
+ | |

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

</haskell> |
</haskell> |
||

− | Well, a natural number ''p'' is a prime number iff it is larger than '''1''' and no natural number ''n'' with ''n >= 2'' and ''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 remainder concerning the division of ''p'' by ''n''. |
+ | 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 natural ''primes'' <= sqrt P. |
+ | However, we don't actually need to check all natural numbers ''<= sqrt P''. We need only check the ''primes <= sqrt P'': |

<haskell> |
<haskell> |
||

-- Infinite list of all prime numbers |
-- Infinite list of all prime numbers |
||

− | allPrimes :: [Int] |
+ | {-# OPTIONS_GHC -O2 -fno-cse #-} |

− | allPrimes = filter (isPrime) [2..] |
+ | candidateFactors p = let z = floor $ sqrt $ fromIntegral p + 1 in |

+ | takeWhile (<= z) primesTME |
||

− | isPrime :: Int -> Bool |
+ | -- tree-merging Eratosthenes sieve |

− | isPrime p |
+ | primesTME = 2 : gaps 3 (join [[p*p,p*p+2*p..] | p <- primes']) |

− | | p < 2 = error "Number too small" |
+ | where |

− | | p == 2 = True |
+ | primes' = 3 : gaps 5 (join [[p*p,p*p+2*p..] | p <- primes']) |

− | | p > 2 = all (\n -> p `mod` n /= 0) (getPrimes sqrtp) |
+ | join ((x:xs):t) = x : union xs (join (pairs t)) |

− | where getPrimes z = takeWhile (<= z) allPrimes |
+ | pairs ((x:xs):ys:t) = (x : union xs ys) : pairs t |

− | sqrtp = floor . sqrt $ fromIntegral p |
+ | gaps k xs@(x:t) | k==x = gaps (k+2) t |

− | </haskell> |
+ | | True = k : gaps (k+2) xs |

− | Note that the mutual dependency of allPrimes and isPrime would result in an infinite loop if we weren't careful. But since we limit our observation of allPrimes to <= sqrt x, we avoid infinite recursion. |
+ | -- 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 |
||

+ | </haskell> |
||

− | While the mutual dependency is interesting, this second version is not necessarily more efficient than the first. Though we avoid checking all natural numbers <= sqrt P in the isPrime method, we instead check the primality of all natural numbers <= sqrt P in the allPrimes definition. |
+ | The tree-merging Eratosthenes sieve here seems to strike a good balance between efficiency and brevity. More at [[Prime numbers]] haskellwiki page. |

## Revision as of 07:58, 31 May 2011

(**) 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*:

-- Infinite list of all prime numbers {-# OPTIONS_GHC -O2 -fno-cse #-} candidateFactors p = let z = floor $ sqrt $ fromIntegral p + 1 in takeWhile (<= z) primesTME -- tree-merging Eratosthenes sieve 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.