From 93b20967b5f9a62a522c37a7dec49ae6037d89c0 Mon Sep 17 00:00:00 2001
From: Samuel Tardieu <sam@rfc1149.net>
Date: Fri, 26 Dec 2008 20:58:46 +0100
Subject: [PATCH] Add Eratosthene sieve as math.primes.erato

---
 extra/math/primes/erato/erato-docs.factor  | 12 +++++++++++
 extra/math/primes/erato/erato-tests.factor |  3 +++
 extra/math/primes/erato/erato.factor       | 23 ++++++++++++++++++++++
 3 files changed, 38 insertions(+)
 create mode 100644 extra/math/primes/erato/erato-docs.factor
 create mode 100644 extra/math/primes/erato/erato-tests.factor
 create mode 100644 extra/math/primes/erato/erato.factor

diff --git a/extra/math/primes/erato/erato-docs.factor b/extra/math/primes/erato/erato-docs.factor
new file mode 100644
index 0000000000..b12ea45052
--- /dev/null
+++ b/extra/math/primes/erato/erato-docs.factor
@@ -0,0 +1,12 @@
+USING: help.markup help.syntax ;
+IN: math.primes.erato
+
+HELP: sieve
+{ $values { "n" "the greatest odd number to consider" } { "arr" "a bit array" } }
+{ $description "Return a bit array containing a primality bit for every odd number between 3 and " { $snippet "n" } " (inclusive). " { $snippet ">index" } " can be used to retrieve the index of an odd number to be tested." } ;
+
+HELP: >index
+{ $values { "n" "an odd number" } { "i" "the corresponding index" } }
+{ $description "Retrieve the index corresponding to the odd number on the stack." } ;
+
+{ sieve >index } related-words
diff --git a/extra/math/primes/erato/erato-tests.factor b/extra/math/primes/erato/erato-tests.factor
new file mode 100644
index 0000000000..917824c9c1
--- /dev/null
+++ b/extra/math/primes/erato/erato-tests.factor
@@ -0,0 +1,3 @@
+USING: bit-arrays math.primes.erato tools.test ;
+
+[ ?{ t t t f t t f t t f t f f t } ] [ 29 sieve ] unit-test
diff --git a/extra/math/primes/erato/erato.factor b/extra/math/primes/erato/erato.factor
new file mode 100644
index 0000000000..f4409038bb
--- /dev/null
+++ b/extra/math/primes/erato/erato.factor
@@ -0,0 +1,23 @@
+USING: bit-arrays kernel math math.functions math.ranges sequences ;
+IN: math.primes.erato
+
+: >index ( n -- i )
+    3 - 2 /i ; inline
+
+: index> ( i -- n )
+    2 * 3 + ; inline
+
+: mark-multiples ( i arr -- )
+    [ dup index> [ + ] keep ] dip
+    [ length 1 - swap <range> f swap ] keep
+    [ set-nth ] curry with each ;
+
+: maybe-mark-multiples ( i arr -- )
+    2dup nth [ mark-multiples ] [ 2drop ] if ;
+
+: init-sieve ( n -- arr )
+    >index 1 + <bit-array> dup set-bits ;
+
+: sieve ( n -- arr )
+    [ init-sieve ] [ sqrt >index [0,b] ] bi
+    over [ maybe-mark-multiples ] curry each ; foldable