This is an automated email from the git hooks/post-receive script.

ppm-guest pushed a commit to annotated tag v0.36
in repository libmath-prime-util-perl.

commit 0d643508f1a93c4ae433019a59fc120e7ac35140
Author: Dana Jacobsen <d...@acm.org>
Date:   Fri Dec 20 19:55:07 2013 -0800

    New sieve code for finding segment composite inc/masks
---
 Changes |   2 +
 sieve.c | 147 +++++++++++++++++++++++++++++++++++++++++++++++++---------------
 sieve.h |  11 -----
 3 files changed, 116 insertions(+), 44 deletions(-)

diff --git a/Changes b/Changes
index 4299680..785ff21 100644
--- a/Changes
+++ b/Changes
@@ -33,6 +33,8 @@ Revision history for Perl module Math::Prime::Util
       lehmer.c.  Really only affects LMOS and Legendre, but it's significant
       for them.  Thanks to Kim Walisch for questioning my earlier decision.
 
+    - Rewrite sieve composite map.  Sieving is now a little faster.
+
 
 0.35  2013-12-08
 
diff --git a/sieve.c b/sieve.c
index 95c90b9..5b77095 100644
--- a/sieve.c
+++ b/sieve.c
@@ -95,30 +95,116 @@ static const unsigned char presieve13[PRESIEVE_SIZE] =
   0x18,0x89,0x08,0x25,0x44,0x22,0x30,0x14,0xc3,0x88,0x86,0x40,0x1a,
   0x28,0x30,0x85,0x09,0x54,0x60,0x43,0x24,0x92,0x81,0x08,0x04,0x70};
 
-#define FIND_COMPOSITE_POS(i,j) \
-  { \
-    UV dlast = d; \
-    do { \
-      d += dinc; \
-      m += minc; \
-      if (m >= 30) { d++; m -= 30; } \
-    } while ( masktab30[m] == 0 ); \
-    wdinc[i] = d - dlast; \
-    wmask[j] = masktab30[m]; \
-  }
-#define FIND_COMPOSITE_POSITIONS(p) \
+static const unsigned char stepdata[8][8][8] = {
+  { {96,65,34,67,36,69,102,47},
+    {64,47,75,106,46,105,77,44},
+    {32,86,33,87,115,45,114,84},
+    {72,46,85,42,73,127,44,123},
+    {96,35,100,79,33,66,37,78},
+    {64,108,34,109,67,47,65,46},
+    {32,76,117,33,118,74,35,87},
+    {32,127,86,45,84,43,82,121},
+  },
+  { {65,34,67,36,69,102,47,96},
+    {105,77,44,64,47,75,106,46},
+    {33,87,115,45,114,84,32,86},
+    {73,127,44,123,72,46,85,42},
+    {33,66,37,78,96,35,100,79},
+    {65,46,64,108,34,109,67,47},
+    {33,118,74,35,87,32,76,117},
+    {121,32,127,86,45,84,43,82},
+  },
+  { {34,67,36,69,102,47,96,65},
+    {106,46,105,77,44,64,47,75},
+    {114,84,32,86,33,87,115,45},
+    {42,73,127,44,123,72,46,85},
+    {66,37,78,96,35,100,79,33},
+    {34,109,67,47,65,46,64,108},
+    {74,35,87,32,76,117,33,118},
+    {82,121,32,127,86,45,84,43},
+  },
+  { {67,36,69,102,47,96,65,34},
+    {75,106,46,105,77,44,64,47},
+    {115,45,114,84,32,86,33,87},
+    {123,72,46,85,42,73,127,44},
+    {35,100,79,33,66,37,78,96},
+    {67,47,65,46,64,108,34,109},
+    {35,87,32,76,117,33,118,74},
+    {43,82,121,32,127,86,45,84},
+  },
+  { {36,69,102,47,96,65,34,67},
+    {44,64,47,75,106,46,105,77},
+    {84,32,86,33,87,115,45,114},
+    {44,123,72,46,85,42,73,127},
+    {100,79,33,66,37,78,96,35},
+    {108,34,109,67,47,65,46,64},
+    {76,117,33,118,74,35,87,32},
+    {84,43,82,121,32,127,86,45},
+  },
+  { {69,102,47,96,65,34,67,36},
+    {77,44,64,47,75,106,46,105},
+    {45,114,84,32,86,33,87,115},
+    {85,42,73,127,44,123,72,46},
+    {37,78,96,35,100,79,33,66},
+    {109,67,47,65,46,64,108,34},
+    {117,33,118,74,35,87,32,76},
+    {45,84,43,82,121,32,127,86},
+  },
+  { {102,47,96,65,34,67,36,69},
+    {46,105,77,44,64,47,75,106},
+    {86,33,87,115,45,114,84,32},
+    {46,85,42,73,127,44,123,72},
+    {78,96,35,100,79,33,66,37},
+    {46,64,108,34,109,67,47,65},
+    {118,74,35,87,32,76,117,33},
+    {86,45,84,43,82,121,32,127},
+  },
+  { {47,96,65,34,67,36,69,102},
+    {47,75,106,46,105,77,44,64},
+    {87,115,45,114,84,32,86,33},
+    {127,44,123,72,46,85,42,73},
+    {79,33,66,37,78,96,35,100},
+    {47,65,46,64,108,34,109,67},
+    {87,32,76,117,33,118,74,35},
+    {127,86,45,84,43,82,121,32},
+  },
+};
+
+static const int wheelmap[30] =
+  {0,0,0,0,0,0,0,1,0,0,0,2,0,3,0,0,0,4,0,5,0,0,0,6,0,0,0,0,0,7};
+static const int wheel2xmap[30] =     /* (2*p)%30 => 2,14,22,26,4,8,16,28 */
+  {0,0,0,0,4,0,0,0,5,0,0,0,0,0,1,0,6,0,0,0,0,0,2,0,0,0,3,0,7,0};
+/*     2   4       8          14  16          22      26  28     (2*p)%30 */
+
+#define FIND_COMPOSITE_POSITIONS(d, m, p) \
   do { \
-    FIND_COMPOSITE_POS(0,1) \
-    FIND_COMPOSITE_POS(1,2) \
-    FIND_COMPOSITE_POS(2,3) \
-    FIND_COMPOSITE_POS(3,4) \
-    FIND_COMPOSITE_POS(4,5) \
-    FIND_COMPOSITE_POS(5,6) \
-    FIND_COMPOSITE_POS(6,7) \
-    FIND_COMPOSITE_POS(7,0) \
-    d -= p; \
+    int v; \
+    UV dinc = (2*p) / 30; \
+    UV minc = (2*p) - dinc*30; \
+    const unsigned char* steps = stepdata [wheelmap[m]] [wheel2xmap[minc]]; \
+    v = steps[0]; wdinc[0] = dinc*(v>>5)+((v>>3)&0x3); wmask[0] = 1<<(v&0x7); \
+    v = steps[1]; wdinc[1] = dinc*(v>>5)+((v>>3)&0x3); wmask[1] = 1<<(v&0x7); \
+    v = steps[2]; wdinc[2] = dinc*(v>>5)+((v>>3)&0x3); wmask[2] = 1<<(v&0x7); \
+    v = steps[3]; wdinc[3] = dinc*(v>>5)+((v>>3)&0x3); wmask[3] = 1<<(v&0x7); \
+    v = steps[4]; wdinc[4] = dinc*(v>>5)+((v>>3)&0x3); wmask[4] = 1<<(v&0x7); \
+    v = steps[5]; wdinc[5] = dinc*(v>>5)+((v>>3)&0x3); wmask[5] = 1<<(v&0x7); \
+    v = steps[6]; wdinc[6] = dinc*(v>>5)+((v>>3)&0x3); wmask[6] = 1<<(v&0x7); \
+    v = steps[7]; wdinc[7] = dinc*(v>>5)+((v>>3)&0x3); wmask[7] = 1<<(v&0x7); \
   } while (0)
 
+/* 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}
+};
+
+
 static void sieve_prefill(unsigned char* mem, UV startd, UV endd)
 {
   UV nbytes = endd - startd + 1;
@@ -141,8 +227,6 @@ static void sieve_prefill(unsigned char* mem, UV startd, UV 
endd)
   }
 }
 
-
-
 /* Wheel 30 sieve.  Ideas from Terje Mathisen and Quesada / Van Pelt. */
 unsigned char* sieve_erat30(UV end)
 {
@@ -164,15 +248,14 @@ unsigned char* sieve_erat30(UV end)
 
   limit = isqrt(end);  /* prime*prime can overflow */
   for (prime = 17; prime <= limit; prime = next_prime_in_sieve(mem,prime)) {
-    UV d = (prime*prime)/30;
-    UV m = (prime*prime) - d*30;
-    UV dinc = (2*prime)/30;
-    UV minc = (2*prime) - dinc*30;
+    UV p2 = prime*prime;
+    UV d = p2 / 30;
+    UV m = p2 - d*30;
     UV wdinc[8];
     unsigned char wmask[8];
 
     /* Find the positions of the next composites we will mark */
-    FIND_COMPOSITE_POSITIONS(prime);
+    FIND_COMPOSITE_POSITIONS(d, m, prime);
 #if 0
     assert(d == ((prime*prime)/30));
     assert(d < max_buf);
@@ -257,16 +340,14 @@ int sieve_segment(unsigned char* mem, UV startd, UV endd)
         do { p2 += 2*p; } while (masktab30[p2%30] == 0);
       } while ( (p2 <= endp) && (p2 >= startp) );
 #else
-      UV d = (p2)/30;
-      UV m = (p2) - d*30;
-      UV dinc = (2*p)/30;
-      UV minc = (2*p) - dinc*30;
+      UV d = p2 / 30;
+      UV m = p2 - d*30;
       UV wdinc[8];
       unsigned char wmask[8];
       UV offset_endd = endd - startd;
 
       /* Find the positions of the next composites we will mark */
-      FIND_COMPOSITE_POSITIONS(p);
+      FIND_COMPOSITE_POSITIONS(d, m, p);
       d -= startd;
       /* Unrolled inner loop for marking composites */
       while ( (d+p) <= offset_endd ) {
diff --git a/sieve.h b/sieve.h
index 7fd7610..aeb3bdb 100644
--- a/sieve.h
+++ b/sieve.h
@@ -29,17 +29,6 @@ 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) {

-- 
Alioth's /usr/local/bin/git-commit-notice on 
/srv/git.debian.org/git/pkg-perl/packages/libmath-prime-util-perl.git

_______________________________________________
Pkg-perl-cvs-commits mailing list
Pkg-perl-cvs-commits@lists.alioth.debian.org
http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/pkg-perl-cvs-commits

Reply via email to