Skip to content

Instantly share code, notes, and snippets.

@WillNess

WillNess/Turner's to Bird's.hs

Last active December 20, 2015 00:19
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save WillNess/6040471 to your computer and use it in GitHub Desktop.
Save WillNess/6040471 to your computer and use it in GitHub Desktop.
Download ZIP
Turner's to Bird's
Raw
Turner's to Bird's.hs
sieve (p:xs) = p : sieve [x | x <- xs, rem x p /= 0] -- (0) Turner's
ps = sieve [2..]
------------------------------------------------------------------------
ps = 2 : sieve [3..] ps
sieve xs (p:pt) | (h,t) <- span (< p*p) xs = -- (1)
h ++ sieve [x | x<-t, rem x p /= 0] pt -- Postponed Turner's
------------------------------------------------------------------------
-- rem x p /= 0 === not $ ordElem x [p,p+p..]
-- [x | x <- xs, not $ ordElem x [p,p+p..]] === diff xs [p,p+p..]
-- [x | x <- xs, (x /=).head.snd.span(< x) $ [p,p+p..]] ... -- no repeats in the
-- diff xs [p,p+p..] -- ascending xs, ys
ordzip f g h a b = go a b where go a@(x:t) b@(y:r) = case compare x y of
LT -> f x (go t b) ; EQ -> g x (go t r) ; GT -> h y (go a r)
skip a b = b
union a b = ordzip (:) (:) (:) a b
diff a b = ordzip (:) skip skip a b
------------------------------------------------------------------------
sieve (p:xs) = p : sieve (diff xs [p,p+p..]) -- (0b) basic
ps = sieve [2..] -- Eratosthenes
------------------------------------------------------------------------
-- The Sieve of Eratosthenes is *dual* to Trial Division - where the latter
-- searches for a number's divisor among primes preceding its square root, the former
-- generates a prime's multiples starting from its square!
-- TD suspects every number to be a composite and must exhaustively prove it's a prime;
-- SoE assumes each number is a prime and just enumerates its multiples (to be eliminated)
------------------------------------------------------------------------
ps = 2 : sieve [3..] ps
sieve xs (p:pt) | (h,t) <- span (< p*p) xs = -- (1) -> (3)
h ++ sieve (diff t [p*p, p*p+p..]) pt -- Postponed Eratosthenes'
ps = 2 : sieve 3 ps []
sieve x (p:pt) cs | q<-p*p, (h,t) <- span (< q) cs = -- (4)
diff [x..q-1] h ++ sieve q pt (union t [q, q+p..]) -- Segmentwise Postponed
-- Eratosthenes'
ps = 2 : 3 : diff [4..] (sieve ps [])
sieve (p:pt) cs | q<-p*p, (h,t) <- span (< q) cs = -- (5)
h ++ sieve pt (union t [q, q+p..]) -- Combined Composites
ps = 2 : diff [3..] (foldr (\p-> (p*p:).union [p*p+p,p*p+2*p..]) [] ps) -- (6) Bird's
------------------------------------------------------------------------
ps = 2 : sieve 3 (inits ps) ps
sieve x (fs:ft) (p:t) | q <- p*p = -- (2)
[x | x<-[x..q-1], all ((/=0).rem x) fs] ++ sieve q ft t -- Segmentwise generative
-- Turner's
ps = 2 : sieve 3 (inits ps) ps
sieve x (fs:ft) (p:t) | q <- p*p = -- (2),(4) -> (7)
(`diff` foldr (\p-> let s=(x`div`p)*p in union [s,s+p..q-1]) [] fs) -- Segmentwise generative
[x..q-1] ++ sieve q ft t -- Eratosthenes'
(3)->(4) idea from https://gist.github.com/WillNess/5659214 i.e. from http://stackoverflow.com/q/16271592
foldi f (x:t) = f x . foldi f . unfoldr(\(a:b:c)-> Just(f a b,c)) $ t
_U ((x:xs):t) = (x:) . union xs . _U . unfoldr(\(a:b:c)-> Just(union a b,c)) $ t
-- _U = foldi (\(x:xs)-> (x:).union xs)
------------------------------------------------------------------------
_Y g = g (_Y g) -- = g . g . g . g . ...
ps = _Y $ (2 :) . diff [3..] . foldr (\p-> (p*p:).union [p*p+p,p*p+2*p..]) []
= (2:) . _Y $ (3:) . diff [5,7..] . foldi (\(x:xs)->(x:).union xs)
. map (\p-> [p*p, p*p+2*p ..])
= (2:) . _Y $ (3:) . minus (tail . scanl (+) 3 $ cycle [2])
. foldi (\(x:xs)->(x:).union xs)
-- . map (\p-> map (p*) $ scanl (+) p (cycle [2]))
. map (\p-> scanl (\a d->a+p*d) (p*p) (cycle [2]))
-- ((.(p*)).(+))
-- (curry(uncurry (+) . second (p*)))
= ([2,3,5,7]++) . _Y $ (11:) . tail . minus (scanl (+) 11 wh11)
. foldi (\(x:xs)-> (x:) . union xs)
. map (\(w,p)-> scanl (\c d-> c + p*d) (p*p) w)
. equalsBy snd (tails wh11 `zip` scanl (+) 11 wh11)
-- general code: https://ideone.com/nuoLUE foldi: http://www.haskell.org/
-- gaps&hits: https://ideone.com/vkXCXt /haskellwiki/Fold#Tree-like_folds
-- http://stackoverflow.com/questions/20667386/
-- racket-sieve-of-eratosthenes-using-stream/20667487#20667487
primes = sieve [2..]
sieve (p:r) = p : sieve [n | n <- r, rem n p > 0]
-- do less, get done more:
sieves s@(p:r) = s : sieves [n | n <- r, rem n p > 0]
primesTo m = map head h ++ takeWhile (<= m) s where
(h,(s:_))=span ((<= m).(^2).head) $ sieves [2..]
-- ^^^^ the whole difference from the above
---
primes = map head $ iterate (\(p:xs)-> filter ((> 0).(`rem` p)) xs) [2..]
primesTo m = (++).map head.fst <*> takeWhile (<= m).head.snd $ -- (intrp'd, GHCi)
span ((<= m).(^2).head) $ iterate (\(p:xs)-> filter ((> 0).(`rem` p)) xs) [2..] -- 1.44, 3.75,
-- ^1.38
primes = _Y $ concatMap fst . iterate (\(_,(p:ps,xs))-> second ((ps,). -- 10k:1.20, 20k:3.21,
filter ((> 0).(`rem`p))) $ span (< p*p) xs) . (\ps-> ([2],(ps,[3..]))) -- ^1.42
-- [[2],[3],[5,7],[11,13,17,19,23],[29,31,37,41,43,47],[53,59,61,67,71,73,79,83,89,
-- 97,101,103,107,109,113],[127,131,137,139,149,151,157,163,167],[173,179,181,191,...
eratos = 2 : _Y ((3:).diff [5,7..].foldi (\(x:xs)->(x:).union xs) -- 1.90, 4.44
.map (\p->[p*p,p*p+2*p..])) -- (/3.7x : 0.51, 1.20 )
_Y g = g (_Y g) -- ^1.22
foldi f (x:xs) = f x (foldi f (pairs xs)) where pairs(a:b:r)=f a b:pairs r
ordzip a b = g a b where { g a@(x:r) b@(y:q) -- http://ideone.com/1apwzA
| x<y = (x,0):g r b | y<x = (0,y):g a q | otherwise = (x,y):g r q } -- high const.f.impl/exec-spec;
diff xs ys = [x | (x,y)<- ordzip xs ys, x/=0 && y==0 ] -- real one runs 3.7x faster
meet xs ys = [y | (x,y)<- ordzip xs ys, x/=0 && y/=0 ] -- (cmpld): (w/ gaps+joinT)
union xs ys = [z | (x,y)<- ordzip xs ys, let z=max x y] -- http://ideone.com/3b9cNM
--- (x/=0 || y/=0) always holds
-- the shortest (very slow)
otherPrimes = nubBy (((>1).).gcd) [2..]
---
also,
import Control.Applicative (<|>) -- or `mplus`
import Data.Maybe (listToMaybe, fromJust)
___________________________________________________________
ordzip a@(x:t) b@(y:r) = case compare x y of
LT -> (Just x, Nothing) : ordzip t b -- actually, Data.These.This x
GT -> (Nothing, Just y) : ordzip a r -- That y
_ -> (Just x, Just y) : ordzip t r -- EQ -- These x y
___________________________________________________________
ordzip a b = (`unfoldr` (a,b)) (Just . (\(a@(x:t), b@(y:r)) ->
case compare x y of LT -> ( (Just x, Nothing), (t,b) )
GT -> ( (Nothing, Just y), (a,r) )
{- EQ -} _ -> ( (Just x, Just y), (t,r) ) ))
___________________________________________________________
ordzip a@(x:t) b@(y:r) = (u,v) : ordzip (maybe a (const t) u)
(maybe b (const r) v)
where (u,v)=(listToMaybe [x | x <= y], listToMaybe [y | y <= x])
___________________________________________________________
diff xs ys = [x | (Just x,Nothing) <- ordzip xs ys] -- catThis .: ordzip
meet xs ys = [y | (Just _,Just y) <- ordzip xs ys] -- (map snd .) catThese .: ordzip
union xs ys = [z | (u,v) <- ordzip xs ys, let Just z = u <|> v] -- map (these id id const) .: ordzip
-- [fromJust (u <|> v) | (u,v) <- ordzip xs ys]
-- map (fromJust . uncurry (<|>)) $ ordzip xs ys
===========================================================
ordzipBy f a@(x:t) b@(y:r) = (u,v) : ordzipBy f (maybe a (const t) u)
(maybe b (const r) v)
where (u,v) = let fx = f x in
(listToMaybe [x | fx <= y], listToMaybe [y | y <= fx])
__________________________________________________________
meetBy f xs ys = [x | (Just x,Just _) <- ordzipBy f xs ys]
__________________________________________________________
____________________________________________
ordzipBy k f g h n a b = loop a b where
loop a@(x:t) b@(y:r) = case compare (k x) y of
LT -> f x (loop t b)
EQ -> g x (loop t r)
GT -> h y (loop a r)
loop a b = n a b
ordzip = ordzipBy id
union a b = ordzip (:) (:) (:) (++) a b -- join, fuse
minus a b = ordzip (:) skip skip const a b -- diff
equalsBy k a b
= ordzipBy k skip (:) skip (\a b->[]) a b -- meet, intersection
skip a b = b
-- mapAccumL :: (acc -> x -> (acc, y)) -> acc -> [x] -> (acc, [y]) --
diff xs ys = concat . snd $ mapAccumL g ys xs where
g ys x | (a,c@ ~(b:d)) <- span (< x) ys = if null c || b>x then (c, [x]) else (d,[])
meet xs ys = concat . snd $ mapAccumL g ys xs where
g ys x | (a,c@ ~(b:d)) <- span (< x) ys = if null c || b>x then (c, []) else (d,[x])
union xs ys | (a,b) <- mapAccumL g ys xs = concat b++a where
g ys x | (a,c@ ~(b:d)) <- span (< x) ys = (if null c || b>x then c else d, a ++ [x])
-- Another go at it:
_______________________________________________________________________________________
sieve (p:xs) = p : sieve [x | x <- xs, rem x p > 0] -- (0) Turner's
ps = sieve [2..]
_______________________________________________________________________________________
_Y g = g (_Y g) -- = g . g . g . g . g . ........
noMultOf x n = [n | rem n x > 0]
( primes = foldr (\x r-> x : (noMultOf x =<< r)) [] [2..] )
primes = map head $ iterate (\(x:xs) -> noMultOf x =<< xs) [2..] -- same as Turner's
primes = map head $ scanl (\(x:xs) _-> noMultOf x =<< xs) [2..] [0..] -- same
primes = map head $ scanl (\(_:xs) x-> noMultOf x =<< xs) [2..] [2..] -- not
primes = concat . map fst . tail
$ scanl (\(_,xs) x-> let (h,t)=span (< x*x) xs in
(h, noMultOf x =<< t)) ([],[2..]) [2..] -- not yet
primes = concat . map fst . tail
$ scanl (\(_,xs) p-> let (h,t)=span (< p*p) xs in
(h, noMultOf p =<< t)) ([],[2..]) (2:tail primes) -- here it is
primes = _Y $ concat . map fst . tail -- _Y g = g . g . g ...
. scanl (\(_,xs) p-> let (h,t)=span (< p*p) xs in
(h, noMultOf p =<< t)) ([],[2..]) . (2:) . tail
primes = concat . _Y $ map fst
. scanl (\(_,xs) p-> let (h,t)=span (< p*p) xs in
(h, noMultOf p =<< t)) ([2,3],[5,7..]) . tail . concat
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment