[ccan] [PATCH 5/5] eratosthenes: Implementation of the Sieve of Eratosthenes
David Gibson
david at gibson.dropbear.id.au
Fri Oct 3 00:15:01 EST 2014
This adds a new "eratosthenes" module which implements the standard
Sieve of Eratosthenes algorithm for locating small prime numbers.
---
ccan/eratosthenes/LICENSE | 1 +
ccan/eratosthenes/_info | 45 ++++++++++++
ccan/eratosthenes/eratosthenes.c | 116 +++++++++++++++++++++++++++++++
ccan/eratosthenes/eratosthenes.h | 27 +++++++
ccan/eratosthenes/test/run-incremental.c | 36 ++++++++++
ccan/eratosthenes/test/run.c | 57 +++++++++++++++
6 files changed, 282 insertions(+)
create mode 120000 ccan/eratosthenes/LICENSE
create mode 100644 ccan/eratosthenes/_info
create mode 100644 ccan/eratosthenes/eratosthenes.c
create mode 100644 ccan/eratosthenes/eratosthenes.h
create mode 100644 ccan/eratosthenes/test/run-incremental.c
create mode 100644 ccan/eratosthenes/test/run.c
diff --git a/ccan/eratosthenes/LICENSE b/ccan/eratosthenes/LICENSE
new file mode 120000
index 0000000..dc314ec
--- /dev/null
+++ b/ccan/eratosthenes/LICENSE
@@ -0,0 +1 @@
+../../licenses/LGPL-2.1
\ No newline at end of file
diff --git a/ccan/eratosthenes/_info b/ccan/eratosthenes/_info
new file mode 100644
index 0000000..a09e8d1
--- /dev/null
+++ b/ccan/eratosthenes/_info
@@ -0,0 +1,45 @@
+#include "config.h"
+#include <stdio.h>
+#include <string.h>
+
+/**
+ * eratosthenes - Sieve of Eratosthenes
+ *
+ * This code implements Eratosthenes' Sieve for efficiently finding
+ * small prime numbers (in this context anything less than several
+ * billion is "small").
+ *
+ * Example:
+ * #include <ccan/eratosthenes/eratosthenes.h>
+ *
+ * int main(int argc, char *argv[])
+ * {
+ * struct eratosthenes s;
+ * unsigned long p;
+ *
+ * eratosthenes_init(&s);
+ * eratosthenes_sieve(&s, atol(argv[1]));
+ *
+ * while ((p = eratosthenes_nextprime(&s, p)) != 0) {
+ * printf("%ld\n", p);
+ * }
+ *
+ * return 0;
+ * }
+ *
+ * License: LGPL (v2.1 or any later version)
+ * Author: David Gibson <david at gibsond.dropbear.id.au>
+ */
+int main(int argc, char *argv[])
+{
+ /* Expect exactly one argument */
+ if (argc != 2)
+ return 1;
+
+ if (strcmp(argv[1], "depends") == 0) {
+ printf("ccan/bitmap\n");
+ return 0;
+ }
+
+ return 1;
+}
diff --git a/ccan/eratosthenes/eratosthenes.c b/ccan/eratosthenes/eratosthenes.c
new file mode 100644
index 0000000..7c39c3a
--- /dev/null
+++ b/ccan/eratosthenes/eratosthenes.c
@@ -0,0 +1,116 @@
+/* Licensed under LGPLv2+ - see LICENSE file for details */
+#include <ccan/eratosthenes/eratosthenes.h>
+#include <ccan/bitmap/bitmap.h>
+
+#include <assert.h>
+#include <stdlib.h>
+
+#define VAL_TO_BIT(v) (((v) - 3) / 2)
+#define LIMIT_TO_NBITS(l) ((l > 2) ? ((l) - 2) / 2 : 0)
+
+#define BIT_TO_VAL(b) (((b) * 2) + 3)
+
+void eratosthenes_init(struct eratosthenes *s)
+{
+ s->limit = 0;
+ s->b = NULL;
+}
+
+void eratosthenes_reset(struct eratosthenes *s)
+{
+ if (s->b)
+ free(s->b);
+ eratosthenes_init(s);
+}
+
+static void eratosthenes_once(struct eratosthenes *s, unsigned long limit, unsigned long p)
+{
+ unsigned long n = VAL_TO_BIT(3*p);
+ unsigned long obits = LIMIT_TO_NBITS(s->limit);
+
+ if (obits > n) {
+ n = obits + p - 1 - ((obits - n - 1) % p);
+ }
+
+ assert((BIT_TO_VAL(n) % p) == 0);
+ assert((BIT_TO_VAL(n) / p) > 1);
+
+ while (n < LIMIT_TO_NBITS(limit)) {
+ bitmap_clear_bit(s->b, n);
+ n += p;
+ }
+}
+
+static void eratosthenes_sieve_(struct eratosthenes *s, unsigned long limit)
+{
+ unsigned long p = 3;
+
+ while ((p * p) < limit) {
+ unsigned long n;
+
+ eratosthenes_once(s, limit, p);
+
+ n = bitmap_ffs(s->b, VAL_TO_BIT(p) + 1, LIMIT_TO_NBITS(limit));
+
+ /* We should never run out of primes */
+ assert(n < LIMIT_TO_NBITS(limit));
+
+ p = BIT_TO_VAL(n);
+ }
+}
+
+void eratosthenes_sieve(struct eratosthenes *s, unsigned long limit)
+{
+ if ((limit < 3) || (limit <= s->limit))
+ /* Nothing to do */
+ return;
+
+ if (s->limit < 3)
+ s->b = bitmap_alloc1(LIMIT_TO_NBITS(limit));
+ else
+ s->b = bitmap_realloc1(s->b, LIMIT_TO_NBITS(s->limit),
+ LIMIT_TO_NBITS(limit));
+
+ if (!s->b)
+ abort();
+
+ eratosthenes_sieve_(s, limit);
+
+ s->limit = limit;
+}
+
+bool eratosthenes_isprime(const struct eratosthenes *s, unsigned long n)
+{
+ assert(n < s->limit);
+
+ if ((n % 2) == 0)
+ return (n == 2);
+
+ if (n < 3) {
+ assert(n == 1);
+ return false;
+ }
+
+ return bitmap_test_bit(s->b, VAL_TO_BIT(n));
+}
+
+unsigned long eratosthenes_nextprime(const struct eratosthenes *s, unsigned long n)
+{
+ unsigned long i;
+
+ if ((n + 1) >= s->limit)
+ return 0;
+
+ if (n < 2)
+ return 2;
+
+ if (n == 2)
+ return 3;
+
+ i = bitmap_ffs(s->b, VAL_TO_BIT(n) + 1, LIMIT_TO_NBITS(s->limit));
+ if (i == LIMIT_TO_NBITS(s->limit))
+ /* Reached the end of the sieve */
+ return 0;
+
+ return BIT_TO_VAL(i);
+}
diff --git a/ccan/eratosthenes/eratosthenes.h b/ccan/eratosthenes/eratosthenes.h
new file mode 100644
index 0000000..0ec4692
--- /dev/null
+++ b/ccan/eratosthenes/eratosthenes.h
@@ -0,0 +1,27 @@
+/* Licensed under LGPLv2+ - see LICENSE file for details */
+#ifndef CCAN_ERATOSTHENES_H_
+#define CCAN_ERATOSTHENES_H_
+
+#include "config.h"
+
+#include <stdbool.h>
+
+#include <ccan/bitmap/bitmap.h>
+
+struct eratosthenes {
+ unsigned long limit;
+ bitmap *b;
+};
+
+void eratosthenes_init(struct eratosthenes *s);
+
+void eratosthenes_reset(struct eratosthenes *s);
+
+void eratosthenes_sieve(struct eratosthenes *s, unsigned long limit);
+
+bool eratosthenes_isprime(const struct eratosthenes *s, unsigned long n);
+
+unsigned long eratosthenes_nextprime(const struct eratosthenes *s,
+ unsigned long n);
+
+#endif /* CCAN_ERATOSTHENES_H_ */
diff --git a/ccan/eratosthenes/test/run-incremental.c b/ccan/eratosthenes/test/run-incremental.c
new file mode 100644
index 0000000..ea6b2b9
--- /dev/null
+++ b/ccan/eratosthenes/test/run-incremental.c
@@ -0,0 +1,36 @@
+#include <ccan/eratosthenes/eratosthenes.h>
+#include <ccan/tap/tap.h>
+
+#include <ccan/eratosthenes/eratosthenes.c>
+
+#define LIMIT 500
+
+#define ok_eq(a, b) \
+ ok((a) == (b), "%s [%u] == %s [%u]", \
+ #a, (unsigned)(a), #b, (unsigned)(b))
+
+int main(void)
+{
+ struct eratosthenes s1, s2;
+ unsigned long n;
+
+ /* This is how many tests you plan to run */
+ plan_tests(LIMIT);
+
+ eratosthenes_init(&s1);
+ eratosthenes_sieve(&s1, LIMIT);
+
+ eratosthenes_init(&s2);
+ for (n = 1; n <= LIMIT; n++)
+ eratosthenes_sieve(&s2, n);
+
+ for (n = 0; n < LIMIT; n++)
+ ok1(eratosthenes_isprime(&s1, n)
+ == eratosthenes_isprime(&s2, n));
+
+ eratosthenes_reset(&s1);
+ eratosthenes_reset(&s2);
+
+ /* This exits depending on whether all tests passed */
+ return exit_status();
+}
diff --git a/ccan/eratosthenes/test/run.c b/ccan/eratosthenes/test/run.c
new file mode 100644
index 0000000..a40b32b
--- /dev/null
+++ b/ccan/eratosthenes/test/run.c
@@ -0,0 +1,57 @@
+#include <ccan/eratosthenes/eratosthenes.h>
+#include <ccan/tap/tap.h>
+
+#include <ccan/eratosthenes/eratosthenes.c>
+
+#define LIMIT 500
+
+#define ok_eq(a, b) \
+ ok((a) == (b), "%s [%u] == %s [%u]", \
+ #a, (unsigned)(a), #b, (unsigned)(b))
+
+static bool test_isprime(unsigned long n)
+{
+ int i;
+
+ if (n < 2)
+ return false;
+
+ for (i = 2; i < n; i++)
+ if ((n % i) == 0)
+ return false;
+
+ return true;
+}
+
+static unsigned long test_nextprime(struct eratosthenes *s, unsigned long n)
+{
+ unsigned long i = n + 1;
+
+ while ((i < LIMIT) && !eratosthenes_isprime(s, i))
+ i++;
+
+ return (i >= LIMIT) ? 0 : i;
+}
+
+int main(void)
+{
+ struct eratosthenes s;
+ unsigned long n;
+
+ /* This is how many tests you plan to run */
+ plan_tests(2 * LIMIT);
+
+ eratosthenes_init(&s);
+
+ eratosthenes_sieve(&s, LIMIT);
+
+ for (n = 0; n < LIMIT; n++) {
+ ok_eq(eratosthenes_isprime(&s, n), test_isprime(n));
+ ok_eq(eratosthenes_nextprime(&s, n), test_nextprime(&s, n));
+ }
+
+ eratosthenes_reset(&s);
+
+ /* This exits depending on whether all tests passed */
+ return exit_status();
+}
--
1.9.3
More information about the ccan
mailing list