# Prime numbers

### From HaskellWiki

In mathematics, a *prime number* (or a *prime*) is a natural number which has exactly two distinct natural number divisors: 1 and itself. The smallest prime is thus 2.

## Contents |

## 1 Prime Number Resources

- At Wikipedia:

- HackageDB packages:
- Numbers: An assortment of number theoretic functions.
- NumberSieves: Number Theoretic Sieves: primes, factorization, and Euler's Totient.
- primes: Efficient, purely functional generation of prime numbers.

## 2 Finding Primes

Any natural number is representable as a product of powers of its prime factors, and so a prime has no prime divisors other than itself. That means that starting with 2, *for each* newly found *prime* we can *eliminate* from the rest of the numbers *all such numbers* that have this prime as their divisor, giving us the *next available* number as next prime. This is known as *sieving* the natural numbers, so that in the end what we are left with are just primes.

### 2.1 Filtering To Keep the Prime Numbers In

#### 2.1.1 The Classic Turner's Sieve

Attributed to David Turner *(SASL Language Manual, 1975)*, the following is a direct translation of that idea, generating a list of all prime numbers:

primes :: [Integer] primes = sieve [2..] where sieve (p:xs) = p : sieve [x | x <- xs, x `mod` p /= 0] -- or: filter ((/=0).(`mod`p)) xs

This should only be considered a *specification*, not a *code*. When run as is, it is *extremely inefficient* because it starts up the filters prematurely, immediately after each prime, instead of only after the prime's square has been reached. To be admitted as prime, *each number* will be *tested for divisibility* here by all its preceding primes, while just those not greater than its square root would suffice. This means that e.g. to find the **1001**st prime (`7927`

), **1000** filters are used, when in fact just the first **24** are needed (up to `89`

's filter only).

So this in effect creates a cascade of nested filters in front of the infinite numbers supply, and in *extremely premature* fashion at that. One way of fixing that would be to *postpone* the creation of filters until the right moment, by decoupling the primes supply from the numbers supply.

#### 2.1.2 Postponed Filters Sieve

Here each filter's creation is postponed until the very moment it's needed.

primes :: [Integer] primes = 2: 3: sieve (tail primes) [5,7..] where sieve (p:ps) xs = h ++ sieve ps [x | x <- t, x `rem` p /= 0] -- or: filter ((/=0).(`rem`p)) t where (h,~(_:t)) = span (< p*p) xs

This can be seen as essential framework for all the code to come. It only tests odd numbers, and only by the primes that are needed, for each *numbers span* between successive squares of primes. To find the **1001**st prime, the divisibility test is performed by only **24** nested filters corresponding to the first **24** odd primes.

Whereas Turner's sieve exhibits near O(*n*^{2}) behavior (*n* referring to the number of primes produced), this one exhibits near O(*n*^{1.5}) behavior, with an orders-of-magnitude speedup.

#### 2.1.3 Odd numbers, by Trial Division

This is also good for generating a few 100,000s primes (when GHC-compiled as well):

primes :: [Integer] primes = 2: 3: filter isPrime [5,7..] where isPrime n = all (notDivs n) $ takeWhile (\p-> p*p <= n) (tail primes) notDivs n p = n `mod` p /= 0

Instead of relying on nested filters, it tests each odd number by an explicit list of all the needed prime factors. But for each number tested it re-fetches this list *anew* which will be *the same* for the increasing spans of numbers between the successive squares of primes.

#### 2.1.4 Generated Spans, by Nested Filters

The other way to go instead of concentrating on the numbers supply, is to directly work on the successive spans between the primes squares.

This version is a bit faster still, creating *158,000 primes* (again, GHC-compiled) in the same time as the postponed filters does 100,000 primes:

primes :: [Integer] primes = 2: 3: sieve [] (tail primes) 3 where notDivsBy d n = n `mod` d /= 0 sieve ds (p:ps) x = foldr (filter . notDivsBy) [x+2,x+4..p*p-2] ds ++ sieve (p:ds) ps (p*p)

This one explicitly maintains the list of primes needed for testing each odds span between successive primes squares, which it explicitly generates. But it tests with nested `filter`

s, which it repeatedly recreates.

#### 2.1.5 Generated Spans, by List of Primes

The list of primes needed to test each range of odds is actually just the prefix of the primes list itself, of known length, and need not be specifically generated at all. Combined with one-call testing by the explicit list of primes, and direct generation of odds between the successive primes squares, this calculates only the bare essentials and makes no extra recalculations:

primes :: [Integer] primes = 2: 3: sieve 0 (tail primes) 3 sieve k (p:ps) x = [n | n <- [x+2,x+4..p*p-2], and [n`rem`p/=0 | p <- fs]] -- or: all ((>0).(n`rem`)) fs ++ sieve (k+1) ps (p*p) where fs = take k (tail primes)

It produces about *222,000 primes* in the same amount of time, and is good for creating about a million first primes, compiled.

The reason to have `sieve`

function available separately too is that it can also be used to produce primes above a given number, as in

primesFrom m = sieve (length h) ps $ m`div`2*2+1 where (h,(_:ps)) = span (<= (floor.sqrt.fromIntegral) m) primes

It can thus produce a few primes e.g. above `239812076741689`

, which is a square of the millionth odd prime, without having to compute all the preceding primes (which would number in trillions).

### 2.2 Getting the Composite Numbers Out

The divisibility testing too should be considered a specification (as in *no multiples of p*), and not a code per se, because although testing *composites* is cheap (as most of them will have small factors, so the test is soon aborted), testing *prime* numbers is costly, and is to be avoided.

#### 2.2.1 Euler's Sieve

With each found prime Euler's sieve removes all its multiples *in advance* so that at each step the list to process is guaranteed to have *no multiples* of any of the preceding primes in it (consists only of numbers *coprime* with all the preceding primes) and thus starts with the next prime:

import Data.List.Ordered (minus) euler xs@(p:xt) = p : euler (xt `minus` map (p*) xs) primes = euler [2..]

This code is extremely inefficient, running at about O(*n*^{2.4}) complexity, and should be regarded a *specification* only. It works very hard trying to avoid double hits on multiples, which are relatively few and cheap to deal with in the first place. Turner's sieve could also be seen as *its rendition*, substituting the built-in `filter`

for `minus`

- computationally different but producing the same results.

But the situation can be improved using a different list representation, for generating lists not from a last element and an increment, but rather a last span and an increment, which entails a set of helpful equivalences:

{- fromElt (x,i) = x : fromElt (x + i,i) === iterate (+ i) x [n..] === fromElt (n,1) -} fromSpan (xs,i) = xs ++ fromSpan (map (+ i) xs,i) {- === concat $ iterate (map (+ i)) xs [n..] === fromSpan ([n],1) fromSpan (p:xt,i) === p : fromSpan (xt ++ [p + i], i) fromSpan (xs,i) `minus` fromSpan (ys,i) === fromSpan (xs `minus` ys, i) map (p*) (fromSpan (xs,i)) === fromSpan (map (p*) xs, p*i) fromSpan (xs,i) === forall (p > 0). fromSpan (concat $ take p $ iterate (map (+ i)) xs, p*i) -} eulerS () = iterate eulerStep ([2],1) eulerStep (xs@(p:_), i) = ( (tail . concat . take p . iterate (map (+ i))) xs `minus` map (p*) xs, p*i ) {- > mapM_ print $ take 4 $ eulerS () ([2],1) ([3],2) ([5,7],6) ([7,11,13,17,19,23,29,31],30) -}

This is where the notion of *wheel* comes from, naturally and elegantly. For each wheel specification `w@((p:_),_)`

produced by `eulerStep`

, the numbers in `(fromSpan w)`

up to *p*^{2} are all primes too, so that

eulerPrimesTo m = if m > 1 then go ([2],1) else [] where go w@((p:_), _) | m < p*p = takeWhile (<= m) (fromSpan w) | True = p : go (eulerStep w)

This runs at about O(*n*^{1.5}) complexity, for `n`

primes produced.

#### 2.2.2 Postponed Multiples Removal

All the filtering versions thus far try to *keep the primes* among all numbers by testing *each number* in isolation. Instead, we can just *remove* the prime's *multiples* in advance. We gain in speed because we now get the primes *for free*, after all the multiples are removed on a particular span, *without* performing any divisibility tests *at all*, with this simple variation of the Postponed Filters sieve:

primes :: [Integer] primes = 2: 3: sieve (tail primes) [5,7..] where sieve (p:ps) xs = h ++ sieve ps (t `minus` tail [q,q+2*p..]) where (h,~(_:t)) = span (< q) xs q = p*p

This is similar to Euler's sieve, in that it removes multiples of each prime in advance, but will generate some of the multiples more than once, like e.g. `3*15`

and `5*9`

, unlike the true Euler's sieve which won't generate the latter.

##### 2.2.2.1 Merged Multiples Removal Sieve

The left-associating `(...((xs-a)-b)-...)`

above buries the *initial numbers supply* `xs`

(here, `[5,7..]`

list) deeper and deeper on the left, with time. But it is the same as `(xs-(a+b+...))`

, and so we can just remove from it the infinite primes multiples lists (all united into one without duplicates), each starting at its seeding prime's square (arriving at what is essentially equivalent to the Richard Bird's code from Melissa O'Neill's article). This way we needn't explicitly jump to a prime's square because it's where its multiples start anyway:

import Data.List.Ordered (minus,union) primes :: [Integer] primes = 2: primes' where primes' = 3: [5,7..] `minus` foldr union' [] [(q,tail [q,q+2*p..]) | p<-primes', let q=p*p] union' (q,qs) xs = q : union qs xs

This code is yet faster. Its main deficiency still is that it creates a *linear* nested merging structure, instead of a *tree-like structure*. Each multiple produced by a prime has to percolate to the top eventually, so it's better to make its path shorter. It'll have to go through fewer merge nodes this way.

The *linearity* of structure is imposed by *type asymmetry* of our `(union' :: a -> b -> b)`

function, forcing us into the `1+(1+(1+(1+...)))`

pattern, `+`

standing for `union'`

(which was defined that way to prevent the run-ahead when folding over the infinite list of lists of multiples).

We need to turn it into an associative operation of uniform type `(:: a -> a -> a)`

to be able to freely rearrange the combinations into arbitrary tree-like patterns, as in e.g. `((1+1)+(2+2))+(...)`

etc. The type uniformity is what makes compositionality possible.

That is what the Implicit Heap section below improves on, which is essentially equivalent to using a tree-producing `tfold`

instead of a standard linear `foldr`

, in the following:

##### 2.2.2.2 Treefold Merged Multiples Removal

primes :: [Integer] primes = 2: primes' where primes' = 3:5:7:11:13:17: [19,21..] `minus` (fst.tfold unionSP) [([q],tail [q,q+2*p..]) | p<-primes', let q=p*p] unionSP (a:as,bs) v = (a:x,y) where (x,y)=unionSP (as,bs) v unionSP u@([],b:bs) v@(c:cs,ds) = case compare b c of LT -> (b:x,y) where (x,y)=unionSP ([],bs) v EQ -> (b:x,y) where (x,y)=unionSP ([],bs) (cs,ds) GT -> (c:x,y) where (x,y)=unionSP u (cs,ds) unionSP ([],bs) ([],ds) = ([] ,union bs ds) tfold f xs = go (pairs xs) where go (a:b:c:ys) = f a (f b c) `f` go (pairs ys) pairs (a:b:ys) = f a b : pairs ys

The fold used here creates a `(2+(2+2))+( (4+(4+4))+( (8+(8+8))+( ... )))`

structure, better adjusted for primes multiples production than `1+(2+(4+(8+...)))`

, used by the "Implicit Heap" code, giving it additional 10% speedup.

` unionSP `

is an associative operation, preserving of the invariant such that for a list of multiples `[(a,b),(c,d), ...]`

, it's always the case that ` last a < head b && last a < head c`

. These "split pairs" represent ordered list as a pair of its known and (currently) finite prefix, and the rest of it. Such pairs under `unionSP`

operation form a *monoid*, and we could declare a

`unionSP`

its `tfold`

its This code exhibits approximately O(*n*^{1.20})..O(*n*^{1.15}) local asymptotic behavior (tested interpreted, in GHCi, for 10,000 to 300,000 primes produced). When compared with Melissa O'Neill's PQ code from the ZIP package which was modified to work on odds only as well, it is 3.2x times faster, with used memory reported about 2.5x times smaller.

It can be further improved with the wheel optimization (as seen in Euler's Sieve and described below in Prime Wheels section) :

##### 2.2.2.3 Treefold Merged Multiples, with Wheel

The odds-only feed is exactly what `fromSpan([3],2)`

produces. Using a larger wheel means all the multiples of the starting few primes, and not just of `2`

, will be automatically excluded in advance. Coupled with the `People`

-style, more "flat" and "local" data type (from Implicit Heap section below) which is conceptually the same as split list pair but provides for more immediate access to its first element, it leads to:

{-# OPTIONS_GHC -O2 -fno-cse #-} data A a = A a (A a) | B [a] -- direct encoding for Split List primes :: [Integer] primes = 2:3:5:7:primes' where primes' = roll 11 wheel `minus` tfold add [A [p*p] (B [p*q | q<-rollFrom p]) | p<-primes_] primes_ = h ++ t `minus` tfold add [A [p*p] (B [p*q | q<-rollFrom p]) | p<-primes_] where (h,t) = splitAt 6 (roll 11) rollFrom n = go wheelNums wheel where m = (n-11) `mod` 210 go (x:xs) (w:ws) | x==m = roll (n+w) ws | True = go xs ws wheelNums = roll 0 wheel roll = scanl (+) wheel = 2:4:2:4:6:2:6:4:2:4:6:6:2:6:4:2:6:4:6:8:4:2:4:2: 4:8:6:4:6:2:4:6:2:6:6:4:2:4:6:2:6:4:2:4:2:10:2:10:wheel add (A x xs) v = A x (add xs v) add u@(B(x:xs)) v@(A y ys) = case compare x y of LT -> A x (add (B xs) v) EQ -> A x (add (B xs) ys) GT -> A y (add u ys) add (B xs) (B ys) = B (union xs ys) minus a@(x:xs) b@(A y ys) = case compare x y of LT -> x : minus xs b EQ -> minus xs ys GT -> minus a ys

The double primes feed is introduced here to save memory as per Melissa O'Neill's code, and is dependent on no expression sharing being performed by a compiler. This code runs at comparable speed on *Ideone.com* as immutable arrays version, yet is near constant and very low in its memory use (low MBs for first few million primes).

#### 2.2.3 Multiples Removal on Generated Spans, or Sieve of Eratosthenes

Postponed multiples removal, with work done segment-wise, on the spans of numbers between the successive primes squares, becomes

import Data.List.Ordered (minus,union) primes :: [Integer] primes = 2: 3: sieve [] (tail primes) 3 where sieve fs (p:ps) x = [x+2,x+4..q-2] `minus` foldl union [] mults ++ sieve fs' ps q where q = p*p mults = [[y+s,y+2*s..q] | (s,y) <- fs] fs' = (2*p,q) : zip (map fst fs) (map last mults)

This modifies the Generated Spans sieve to "mark" the odd composites in a given range (instead of testing all the odds' divisibility) by generating them - just as a person performing the original sieve of Eratosthenes would do, marking one by one the multiples of the relevant primes. These composites are independently generated so some will be generated multiple times.

Interpreted under both GHCi and WinHugs, it runs *faster*, takes *less* memory, and has better asymptotic behavior, its performance approximately the same as in the Merged Multiples Removal sieve. The advantage in working with spans explicitly is that this code is easily amendable to using arrays for the composites marking and removal on each *finite* span.

### 2.3 Using Immutable Arrays

#### 2.3.1 Generating Segments of Primes

The removal of multiples on each segment of odds in the Sieve of Eratosthenes can be done by actually marking them in an array, instead of manipulating the ordered lists:

primes :: [Integer] primes = 2: 3: sieve [] (tail primes) 3 where sieve fs (p:ps) x = [i | i <- [x+2,x+4..q-2], a ! i] ++ sieve fs' ps q where q = p*p mults = [[y+s,y+2*s..q] | (s,y) <- fs] fs' = (2*p,q) : zip (map fst fs) (map last mults) a = accumArray (\ b c -> False) True (x+2,q-2) [(i,()) | ms <- mults, i <- ms]

The above code can be further streamlined and made more than twice faster by working with odds only, represented as their offsets in segment arrays:

primes :: [Integer] primes = 2: 3: sieve [] (tail primes) 3 where sieve fs (p:ps) x = [i*2 + x | (i,e) <- assocs a, e] ++ sieve fs' ps (p*p) where q = (p*p-x)`div`2 fs' = (p,0) : [(s, rem (y-q) s) | (s,y) <- fs] a = accumArray (\ b c -> False) True (1,q-1) [(i,()) | (s,y) <- fs, i <- [y+s,y+s+s..q]]

Apparently, arrays are *fast*. The above is the fastest code of all presented so far. When run on *Ideone.com* it is somewhat faster than wheeled treefold in producing first few million primes, but is very much unacceptable in its memory consumption which grows faster than O(*n*), quickly getting into tens and hundreds of MBs.

#### 2.3.2 Calculating Primes Upto a Given Value

primesToN n = 2: [i | i <- [3,5..n], ar ! i] where ar = f 5 $ accumArray (\ a b -> False) True (3,n) [(i,()) | i <- [9,15..n]] f p a | q > n = a | True = if null x then a' else f (head x) a' where q = p*p a'= a // [(i,False) | i <- [q,q+2*p..n]] x = [i | i <- [p+2,p+4..n], a' ! i]

#### 2.3.3 Calculating Primes in a Given Range

primesFromTo a b = (if a<3 then [2] else []) ++ [i | i <- [o,o+2..b], ar ! i] where o = max (if even a then a+1 else a) 3 r = floor.sqrt.fromInteger $ b+1 ar = accumArray (\a b-> False) True (o,b) [(i,()) | p <- [3,5..r] , let q = p*p s = 2*p (n,x) = quotRem (o - q) s q' = if o <= q then q else q + (n + signum x)*s , i <- [q',q'+s..b] ]

Although using odds instead of primes, the array generation is so fast that it is very much feasible and even preferable for quick generation of some short spans of relatively big primes.

### 2.4 Using Mutable Arrays

Using mutable arrays is the fastest but not the most memory efficient way to calculate prime numbers in Haskell.

#### 2.4.1 Using ST Array

This method implements the Sieve of Eratosthenes, similar to how you might do it in C, modified to work on odds only. It is fast, but about linear in memory consumption, allocating one whole array for the sequence produced.

import Control.Monad import Control.Monad.ST import Data.Array.IArray import Data.Array.MArray import Data.Array.ST import Data.Array.Unboxed primesToNA :: Int -> UArray Int Bool primesToNA n = runSTUArray sieve where sieve = do let m = (n-1)`div`2 a <- newArray (1,m) True :: ST s (STUArray s Int Bool) let sr = floor . (sqrt::Double->Double) . fromIntegral $ n+1 forM_ [1..sr`div`2] $ \i -> do let ii = 2*i*(i+1) -- == ((2*i+1)^2-1)`div`2 si <- readArray a i when si $ forM_ [ii,ii+i+i+1..m] $ \j -> writeArray a j False return a primesToN :: Int -> [Int] primesToN n = 2:[i*2+1 | (i,e) <- assocs . primesToNA $n, e]

#### 2.4.2 Bitwise prime sieve with Template Haskell

Count the number of prime below a given 'n'. Shows fast bitwise arrays, and an example of Template Haskell to defeat your enemies.

{-# OPTIONS -O2 -optc-O -XBangPatterns #-} module Primes (nthPrime) where import Control.Monad.ST import Data.Array.ST import Data.Array.Base import System import Control.Monad import Data.Bits nthPrime :: Int -> Int nthPrime n = runST (sieve n) sieve n = do a <- newArray (3,n) True :: ST s (STUArray s Int Bool) let cutoff = truncate (sqrt $ fromIntegral n) + 1 go a n cutoff 3 1 go !a !m cutoff !n !c | n >= m = return c | otherwise = do e <- unsafeRead a n if e then if n < cutoff then let loop !j | j < m = do x <- unsafeRead a j when x $ unsafeWrite a j False loop (j+n) | otherwise = go a m cutoff (n+2) (c+1) in loop ( if n < 46340 then n * n else n `shiftL` 1) else go a m cutoff (n+2) (c+1) else go a m cutoff (n+2) c

And place in a module:

{-# OPTIONS -fth #-} import Primes main = print $( let x = nthPrime 10000000 in [| x |] )

Run as:

$ ghc --make -o primes Main.hs $ time ./primes 664579 ./primes 0.00s user 0.01s system 228% cpu 0.003 total

### 2.5 Implicit Heap

The following is a more efficient prime generator, implementing the sieve of Eratosthenes.

See also the message threads Re: "no-coding" functional data structures via lazyness for more about how merging ordered lists amounts to creating an implicit heap and Re: Code and Perf. Data for Prime Finders for an explanation of thedata People a = VIP a (People a) | Crowd [a] mergeP :: Ord a => People a -> People a -> People a mergeP (VIP x xt) ys = VIP x $ mergeP xt ys mergeP (Crowd xs) (Crowd ys) = Crowd $ merge xs ys mergeP xs@(Crowd (x:xt)) ys@(VIP y yt) = case compare x y of LT -> VIP x $ mergeP (Crowd xt) ys EQ -> VIP x $ mergeP (Crowd xt) yt GT -> VIP y $ mergeP xs yt merge :: Ord a => [a] -> [a] -> [a] merge xs@(x:xt) ys@(y:yt) = case compare x y of LT -> x : merge xt ys EQ -> x : merge xt yt GT -> y : merge xs yt diff xs@(x:xt) ys@(y:yt) = case compare x y of LT -> x : diff xt ys EQ -> diff xt yt GT -> diff xs yt foldTree :: (a -> a -> a) -> [a] -> a foldTree f ~(x:xs) = x `f` foldTree f (pairs xs) where pairs ~(x: ~(y:ys)) = f x y : pairs ys primes, nonprimes :: [Integer] primes = 2:3:diff [5,7..] nonprimes nonprimes = serve . foldTree mergeP . map multiples $ tail primes where multiples p = vip [p*p,p*p+2*p..] vip (x:xs) = VIP x $ Crowd xs serve (VIP x xs) = x:serve xs serve (Crowd xs) = xs

### 2.6 Prime Wheels

The idea of only testing odd numbers can be extended further. For instance, it is a useful fact that every prime number other than 2 and 3 must be of the form 6*k* + 1 or 6*k* + 5. Thus, we only need to test these numbers:

primes :: [Integer] primes = 2:3:primes' where 1:p:candidates = [6*k+r | k <- [0..], r <- [1,5]] primes' = p : filter isPrime candidates isPrime n = all (not . divides n) $ takeWhile (\p -> p*p <= n) primes' divides n p = n `mod` p == 0

Such a scheme to generate candidate numbers first that avoid a given set of primes as divisors is called a **prime wheel**. Imagine that you had a wheel of circumference 6 to be rolled along the number line. With spikes positioned 1 and 5 units around the circumference, rolling the wheel will prick holes exactly in those positions on the line whose numbers are not divisible by 2 and 3.

A wheel can be represented by its circumference and the spiked positions.

data Wheel = Wheel Integer [Integer]

We prick out numbers by rolling the wheel.

roll (Wheel n rs) = [n*k+r | k <- [0..], r <- rs]

The smallest wheel is the unit wheel with one spike, it will prick out every number.

w0 = Wheel 1 [1]

nextSize (Wheel n rs) p = Wheel (p*n) [r' | k <- [0..(p-1)], r <- rs, let r' = n*k+r, r' `mod` p /= 0]

mkWheel ds = foldl nextSize w0 ds

Now, we can generate prime numbers with a wheel that for instance avoids all multiples of 2, 3, 5 and 7.

primes :: [Integer] primes = small ++ large where 1:p:candidates = roll $ mkWheel small small = [2,3,5,7] large = p : filter isPrime candidates isPrime n = all (not . divides n) $ takeWhile (\p -> p*p <= n) large divides n p = n `mod` p == 0

It's a pretty big wheel with a circumference of 210 and allows us to calculate the first 10000 primes in convenient time.

A fixed size wheel is fine, but how about adapting the wheel size while generating prime numbers? See Euler's Sieve above, or the functional pearl titled Lazy wheel sieves and spirals of primes for more.

### 2.7 Using IntSet for a traditional sieve

module Sieve where import qualified Data.IntSet as I -- findNext - finds the next member of an IntSet. findNext c is | I.member c is = c | c > I.findMax is = error "Ooops. No next number in set." | otherwise = findNext (c+1) is -- mark - delete all multiples of n from n*n to the end of the set mark n is = is I.\\ (I.fromAscList (takeWhile (<=end) (map (n*) [n..]))) where end = I.findMax is -- primes - gives all primes up to n primes n = worker 2 (I.fromAscList [2..n]) where worker x is | (x*x) > n = is | otherwise = worker (findNext (x+1) is) (mark x is)

## 3 Testing Primality

See Testing primality.

## 4 External links

- A collection of prime generators; the file "ONeillPrimes.hs" contains one of the fastest pure-Haskell prime generators. WARNING: Don't use the priority queue from that file for your projects: it's broken and works for primes only by a lucky chance.