[libmath-prime-util-perl] 35/72: Tweak sieve to make more efficient for large bases
Partha P. Mukherjee
ppm-guest at moszumanska.debian.org
Thu May 21 18:49:39 UTC 2015
This is an automated email from the git hooks/post-receive script.
ppm-guest pushed a commit to annotated tag v0.32
in repository libmath-prime-util-perl.
commit 3508783526f9810a13d3291d0460f1904196c30d
Author: Dana Jacobsen <dana at acm.org>
Date: Thu Sep 19 08:55:46 2013 -0700
Tweak sieve to make more efficient for large bases
---
sieve.c | 29 +++++++++++++++--------------
sieve.h | 37 ++++++++++++++++++++++++++++++++-----
2 files changed, 47 insertions(+), 19 deletions(-)
diff --git a/sieve.c b/sieve.c
index 7a8de05..89caf16 100644
--- a/sieve.c
+++ b/sieve.c
@@ -178,10 +178,20 @@ unsigned char* sieve_erat30(UV end)
assert(d < max_buf);
assert(prime = (wdinc[0]+wdinc[1]+wdinc[2]+wdinc[3]+wdinc[4]+wdinc[5]+wdinc[6]+wdinc[7]));
#endif
- /* Regular code to mark composites.
- * i = 0;
- * do { mem[d] |= wmask[i]; d += wdinc[i]; i = (i+1)&7; } while (d < max_buf);
+ /* Regular code to mark composites:
+ * i = 0;
+ * do {mem[d] |= wmask[i]; d += wdinc[i]; i = (i+1)&7;} while (d < max_buf);
* Unrolled version: */
+ while ( (d+prime) < max_buf ) {
+ mem[d] |= wmask[0]; d += wdinc[0];
+ mem[d] |= wmask[1]; d += wdinc[1];
+ mem[d] |= wmask[2]; d += wdinc[2];
+ mem[d] |= wmask[3]; d += wdinc[3];
+ mem[d] |= wmask[4]; d += wdinc[4];
+ mem[d] |= wmask[5]; d += wdinc[5];
+ mem[d] |= wmask[6]; d += wdinc[6];
+ mem[d] |= wmask[7]; d += wdinc[7];
+ }
while (1) {
mem[d] |= wmask[0]; d += wdinc[0]; if (d >= max_buf) break;
mem[d] |= wmask[1]; d += wdinc[1]; if (d >= max_buf) break;
@@ -193,7 +203,6 @@ unsigned char* sieve_erat30(UV end)
mem[d] |= wmask[7]; d += wdinc[7]; if (d >= max_buf) break;
}
}
-
return mem;
}
@@ -233,7 +242,8 @@ int sieve_segment(unsigned char* mem, UV startd, UV endd)
if (p2 < startp) p2 += p;
}
/* Bump to next multiple that isn't divisible by 2, 3, or 5 */
- while (masktab30[p2%30] == 0) { p2 += p; }
+ /* while (masktab30[p2%30] == 0) { p2 += p; } */
+ p2 += p * primestepadvance30[m_>>1][p2%30];
/* It is possible we've overflowed p2, so check for that */
if ( (p2 <= endp) && (p2 >= startp) ) {
/* Sieve from startd to endd starting at p2, stepping p */
@@ -255,14 +265,6 @@ int sieve_segment(unsigned char* mem, UV startd, UV endd)
/* Find the positions of the next composites we will mark */
FIND_COMPOSITE_POSITIONS(p);
d -= startd;
-#if 0
- i = 0; /* Mark the composites */
- do {
- mem[d] |= wmask[i];
- d += wdinc[i];
- i = (i+1) & 7;
- } while (d <= offset_endd);
-#else
/* Unrolled inner loop for marking composites */
while ( (d+p) <= offset_endd ) {
mem[d] |= wmask[0]; d += wdinc[0];
@@ -285,7 +287,6 @@ int sieve_segment(unsigned char* mem, UV startd, UV endd)
mem[d] |= wmask[7]; d += wdinc[7]; if (d > offset_endd) break;
}
#endif
-#endif
}
}
END_DO_FOR_EACH_SIEVE_PRIME;
diff --git a/sieve.h b/sieve.h
index 1bbc0b9..ecf0093 100644
--- a/sieve.h
+++ b/sieve.h
@@ -29,6 +29,17 @@ static const unsigned char distancewheel30[30] =
/* Once on the wheel, add this to get to next spot. In p space, not m. */
static const unsigned char wheeladvance30[30] =
{0,6,0,0,0,0,0,4,0,0,0,2,0,4,0,0,0,2,0,4,0,0,0,6,0,0,0,0,0,2};
+/* For advancing in prime steps in the segment sieve */
+static const unsigned char primestepadvance30[15][30] = {
+ {1,0,5,4,3,2,1,0,3,2,1,0,1,0,3,2,1,0,1,0,3,2,1,0,5,4,3,2,1,0}, {0}, {0},
+ {1,0,3,2,1,2,1,0,3,2,1,0,1,0,5,2,1,0,5,0,3,4,1,0,1,4,3,2,3,0}, {0},
+ {1,0,1,4,3,4,1,0,1,2,3,0,1,0,3,2,3,0,1,0,1,2,5,0,5,2,1,2,3,0},
+ {1,0,3,2,1,2,1,0,3,4,1,0,5,0,3,2,1,0,1,0,3,2,3,0,1,4,5,2,1,0}, {0},
+ {1,0,1,2,5,4,1,0,3,2,3,0,1,0,1,2,3,0,5,0,1,4,3,0,1,2,1,2,3,0},
+ {1,0,3,2,1,2,5,0,5,2,1,0,1,0,3,2,3,0,1,0,3,2,1,0,1,4,3,4,1,0}, {0},
+ {1,0,3,2,3,4,1,0,1,4,3,0,5,0,1,2,5,0,1,0,1,2,3,0,1,2,1,2,3,0}, {0},{0},
+ {1,0,1,2,3,4,5,0,1,2,3,0,1,0,1,2,3,0,1,0,1,2,3,0,1,2,3,4,5,0}
+};
#ifdef FUNC_is_prime_in_sieve
static int is_prime_in_sieve(const unsigned char* sieve, UV p) {
@@ -71,17 +82,27 @@ static UV prev_prime_in_sieve(const unsigned char* sieve, UV p) {
/* Useful macros for the wheel-30 sieve array */
#define START_DO_FOR_EACH_SIEVE_PRIME(sieve, a, b) \
{ \
+ const unsigned char* sieve_ = sieve; \
UV p = a; \
UV l_ = b; \
UV d_ = p/30; \
UV m_ = p-d_*30; \
+ UV lastd_ = l_/30; \
m_ += distancewheel30[m_]; \
+ while (d_ <= lastd_ && (sieve_[d_] & masktab30[m_])) { \
+ m_ = nextwheel30[m_]; if (m_ == 1) { d_++; } \
+ } \
p = d_*30 + m_; \
while ( p <= l_ ) { \
- if ((sieve[d_] & masktab30[m_]) == 0)
#define END_DO_FOR_EACH_SIEVE_PRIME \
- m_ = nextwheel30[m_]; if (m_ == 1) { d_++; } \
+ do { \
+ m_ = nextwheel30[m_]; \
+ if (m_ == 1) { \
+ do { d_++; } while (d_ <= lastd_ && sieve_[d_] == 0xFF); \
+ if (d_ > lastd_) break; \
+ } \
+ } while (sieve_[d_] & masktab30[m_]); \
p = d_*30+m_; \
} \
}
@@ -93,6 +114,7 @@ static UV prev_prime_in_sieve(const unsigned char* sieve, UV p) {
UV m_ = 7; \
UV p = a; \
UV l_ = b; \
+ UV lastd_ = l_/30; \
if (get_prime_cache(l_, &sieve_) < l_) { \
release_prime_cache(sieve_); \
croak("Could not generate sieve for %"UVuf, l_); \
@@ -103,10 +125,12 @@ static UV prev_prime_in_sieve(const unsigned char* sieve, UV p) {
d_ = p/30; \
m_ = p-d_*30; \
m_ += distancewheel30[m_]; \
+ while (d_ <= lastd_ && (sieve_[d_] & masktab30[m_])) { \
+ m_ = nextwheel30[m_]; if (m_ == 1) { d_++; } \
+ } \
p = d_*30 + m_; \
} \
- while ( p <= l_ ) { \
- if (p < 7 || !(sieve_[d_] & masktab30[m_]))
+ while ( p <= l_ ) {
#define RETURN_FROM_EACH_PRIME(x) \
do { release_prime_cache(sieve_); return(x); } while (0)
@@ -115,7 +139,10 @@ static UV prev_prime_in_sieve(const unsigned char* sieve, UV p) {
if (p < 7) { \
p += 1 + (p > 2); \
} else { \
- m_ = nextwheel30[m_]; if (m_ == 1) { d_++; } \
+ do { \
+ m_ = nextwheel30[m_]; \
+ if (m_ == 1) { d_++; if (d_ > lastd_) break; } \
+ } while (sieve_[d_] & masktab30[m_]); \
p = d_*30+m_; \
} \
} \
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-perl/packages/libmath-prime-util-perl.git
More information about the Pkg-perl-cvs-commits
mailing list