Skip to content

Instantly share code, notes, and snippets.

@scmu
Created May 12, 2018
Embed
What would you like to do?
Goldbach's conjecture
{-
Richard Bird's quick algorithm for generating prime numbers.
Note that the usual "one-liner" Haskell definition of the
sieve is not an honest implementation of the sieve of Eratosthenes.
The following algorithm by Bird is a more faithful implementation
that is not asymptotically the best, but practically
fast enough. Further more, we can reuse the "union" function.
For more info see: Melissa E. O’Neill, The genuine sieve of
Eratosthenes. Journal of Functional Programming 19(1), 2009.
-}
primes :: [Int]
primes = 2:([3..] `minus` composites)
where composites = unionBy id [multiples p | p <- primes]
multiples n = map (n*) [n..]
(x:xs) `minus` (y:ys) | x < y = x : (xs `minus` (y:ys))
| x == y = xs `minus` ys
| x > y = (x:xs) `minus` ys
-- merging sorted infinite lists of infinite lists.
unionBy :: Ord b => (a -> b) -> [[a]] -> [a]
unionBy f = foldr merge []
where merge (x:xs) ys = x : merge' xs ys
merge' (x:xs) (y:ys) | f x < f y = x : merge' xs (y:ys)
| f x == f y = x : merge' xs ys
| f x > f y = y : merge' (x:xs) ys
-- the main function
goldbach :: Int -> Maybe (Int, Int)
goldbach n = check . head . dropWhile ((< n) . uncurry (+)) $ primepairs ()
where check (p,q) | p + q == n = Just (p,q)
| otherwise = Nothing
-- all pairs of primes, sorted by their sums
primepairs :: () -> [(Int, Int)]
primepairs _ = unionBy (uncurry (+)) (primesMatrix primes)
where primesMatrix (p:ps) = map (\q -> (p, q)) (p:ps) : primesMatrix ps
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment