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

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

commit ede3138f8237707b307356ab9cc3a28685446c95
Author: Dana Jacobsen <d...@acm.org>
Date:   Mon May 27 00:02:20 2013 -0700

    segmented sieve infrastructure, and have forprimes use it.
---
 XS.xs            | 38 +++++++++++++++-------
 sieve.c          | 97 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 sieve.h          |  3 ++
 t/32-iterators.t | 28 +++++++++++++++-
 util.c           |  2 ++
 5 files changed, 155 insertions(+), 13 deletions(-)

diff --git a/XS.xs b/XS.xs
index 3422caa..9bc19c5 100644
--- a/XS.xs
+++ b/XS.xs
@@ -247,6 +247,7 @@ segment_primes(IN UV low, IN UV high);
   PREINIT:
     AV* av = newAV();
   CODE:
+    /* Could rewrite using {start/next/end}_segment_primes functions */
     if ((low <= 2) && (high >= 2)) { av_push(av, newSVuv( 2 )); }
     if ((low <= 3) && (high >= 3)) { av_push(av, newSVuv( 3 )); }
     if ((low <= 5) && (high >= 5)) { av_push(av, newSVuv( 5 )); }
@@ -650,6 +651,9 @@ forprimes (SV* block, IN SV* svbeg, IN SV* svend = 0)
     HV *stash;
     CV *cv;
     SV* svarg;
+    void* ctx;
+    unsigned char* segment;
+    UV seg_base, seg_low, seg_high;
 
     if (!_validate_int(svbeg, 0) || (items >= 3 && !_validate_int(svend,0))) {
       dSP;
@@ -675,29 +679,39 @@ forprimes (SV* block, IN SV* svbeg, IN SV* svend = 0)
       croak("Not a subroutine reference");
     SAVESPTR(GvSV(PL_defgv));
     svarg = newSVuv(0);
+    ctx = start_segment_primes(beg, end, &segment);
     if (!CvISXSUB(cv)) {
       dMULTICALL;
       I32 gimme = G_VOID;
       PUSH_MULTICALL(cv);
-      START_DO_FOR_EACH_PRIME(beg, end) {
-        sv_setuv(svarg, p);
-        GvSV(PL_defgv) = svarg;
-        MULTICALL;
-      } END_DO_FOR_EACH_PRIME
+      while (beg < 7) {
+        beg = (beg <= 2) ? 2 : (beg <= 3) ? 3 : 5; 
+        { sv_setuv(svarg, beg); GvSV(PL_defgv) = svarg; MULTICALL; }
+        beg += 1 + (beg > 2);
+      }
+      while (next_segment_primes(ctx, &seg_base, &seg_low, &seg_high)) {
+        START_DO_FOR_EACH_SIEVE_PRIME( segment, seg_low - seg_base, seg_high - 
seg_base )
+          { sv_setuv(svarg, seg_base + p); GvSV(PL_defgv) = svarg; MULTICALL; }
+        END_DO_FOR_EACH_SIEVE_PRIME
+      }
 #ifdef PERL_HAS_BAD_MULTICALL_REFCOUNT
       if (CvDEPTH(multicall_cv) > 1)
         SvREFCNT_inc(multicall_cv);
 #endif
       POP_MULTICALL;
     } else {
-      START_DO_FOR_EACH_PRIME(beg, end) {
-        dSP;
-        sv_setuv(svarg, p);
-        GvSV(PL_defgv) = svarg;
-        PUSHMARK(SP);
-        call_sv((SV*)cv, G_VOID|G_DISCARD);
-      } END_DO_FOR_EACH_PRIME
+      while (beg < 7) {
+        beg = (beg <= 2) ? 2 : (beg <= 3) ? 3 : 5; 
+        { dSP; sv_setuv(svarg, beg); GvSV(PL_defgv) = svarg; PUSHMARK(SP); 
call_sv((SV*)cv, G_VOID|G_DISCARD); }
+        beg += 1 + (beg > 2);
+      }
+      while (next_segment_primes(ctx, &seg_base, &seg_low, &seg_high)) {
+        START_DO_FOR_EACH_SIEVE_PRIME( segment, seg_low - seg_base, seg_high - 
seg_base )
+        { dSP; sv_setuv(svarg, seg_base + p); GvSV(PL_defgv) = svarg; 
PUSHMARK(SP); call_sv((SV*)cv, G_VOID|G_DISCARD); }
+        END_DO_FOR_EACH_SIEVE_PRIME
+      }
     }
+    end_segment_primes(ctx);
     SvREFCNT_dec(svarg);
 #endif
     XSRETURN_UNDEF;
diff --git a/sieve.c b/sieve.c
index a8dae1a..5e94d5d 100644
--- a/sieve.c
+++ b/sieve.c
@@ -293,3 +293,100 @@ int sieve_segment(unsigned char* mem, UV startd, UV endd)
   release_prime_cache(sieve);
   return 1;
 }
+
+/**************************************************************************/
+
+typedef struct {
+  UV lod;
+  UV hid;
+  UV low;
+  UV high;
+  UV endp;
+  UV segment_size;
+  unsigned char* segment;
+  unsigned char* base;
+} segment_context_t;
+
+/*
+ * unsigned char* segment;
+ * UV seg_base, seg_low, seg_high;
+ * void* ctx = start_segment_primes(low, high, &segment);
+ * while (beg < 7) {
+ *   beg = (beg <= 2) ? 2 : (beg <= 3) ? 3 : 5;
+ *   .... with beg ....
+ *   beg += 1 + (beg > 2);
+ * }
+ * while (next_segment_primes(ctx, &seg_base, &seg_low, &seg_high)) {
+ *   START_DO_FOR_EACH_SIEVE_PRIME( segment, seg_low - seg_base, seg_high - 
seg_base )
+ *     .... with seg_base + p ....
+ *   END_DO_FOR_EACH_SIEVE_PRIME
+ * }
+ * end_segment_primes(ctx);
+ */
+
+void* start_segment_primes(UV low, UV high, unsigned char** segmentmem)
+{
+  segment_context_t* ctx;
+
+  MPUassert( high >= low, "start_segment_primes bad arguments");
+  New(0, ctx, 1, segment_context_t);
+  ctx->low = low;
+  ctx->high = high;
+  ctx->lod = low / 30;
+  ctx->hid = high / 30;
+  ctx->endp = (ctx->hid >= (UV_MAX/30))  ?  UV_MAX-2  :  30*ctx->hid+29;
+
+  ctx->segment = get_prime_segment( &(ctx->segment_size) );
+  if (ctx->segment == 0)
+    croak("start_segment_primes: Could not get segment");
+  *segmentmem = ctx->segment;
+
+  ctx->base = sieve_erat30( isqrt(ctx->endp)+1 );
+  if (ctx->base == 0)
+    croak("start_segment_primes: Could not get base");
+
+  return (void*) ctx;
+}
+
+int next_segment_primes(void* vctx, UV* base, UV* low, UV* high)
+{
+  UV seghigh_d, range_d;
+  segment_context_t* ctx = (segment_context_t*) vctx;
+
+  if (ctx->lod > ctx->hid) return 0;
+
+  seghigh_d = ((ctx->hid - ctx->lod) < ctx->segment_size)
+            ? ctx->hid
+           : (ctx->lod + ctx->segment_size - 1);
+  range_d = seghigh_d - ctx->lod + 1;
+  *low = ctx->low;
+  *high = (seghigh_d == ctx->hid) ? ctx->high : (seghigh_d*30 + 29);
+  *base = ctx->lod * 30;
+
+  MPUassert( seghigh_d >= ctx->lod, "next_segment_primes: highd < lowd");
+  MPUassert( range_d <= ctx->segment_size, "next_segment_primes: range > 
segment size");
+
+  if (sieve_segment(ctx->segment, ctx->lod, seghigh_d) == 0) {
+    croak("Could not segment sieve from %"UVuf" to %"UVuf, *base+1, *high);
+  }
+
+  ctx->lod += range_d;
+  ctx->low = *high + 2;
+  
+  return 1;
+}
+
+void end_segment_primes(void* vctx)
+{
+  segment_context_t* ctx = (segment_context_t*) vctx;
+  MPUassert(ctx != 0, "end_segment_primes given a null pointer");
+  if (ctx->segment != 0) {
+    release_prime_segment(ctx->segment);
+    ctx->segment = 0;
+  }
+  if (ctx->base != 0) {
+    Safefree(ctx->base);
+    ctx->base = 0;
+  }
+  Safefree(ctx);
+}
diff --git a/sieve.h b/sieve.h
index 4dddde0..1bbc0b9 100644
--- a/sieve.h
+++ b/sieve.h
@@ -6,6 +6,9 @@
 
 extern unsigned char* sieve_erat30(UV end);
 extern int sieve_segment(unsigned char* mem, UV startd, UV endd);
+extern void* start_segment_primes(UV low, UV high, unsigned char** segmentmem);
+extern int next_segment_primes(void* vctx, UV* base, UV* low, UV* high);
+extern void end_segment_primes(void* vctx);
 
 
 static const UV wheel30[] = {1, 7, 11, 13, 17, 19, 23, 29};
diff --git a/t/32-iterators.t b/t/32-iterators.t
index 202ad66..bfc285b 100644
--- a/t/32-iterators.t
+++ b/t/32-iterators.t
@@ -8,7 +8,8 @@ use Math::Prime::Util qw/primes forprimes prime_iterator/;
 my $use64 = Math::Prime::Util::prime_get_config->{'maxbits'} > 32;
 
 plan tests => 8 + 5
-            + 3 + 7;
+            + 3 + 7
+            + 2;
 
 ok(!eval { forprimes { 1 } undef; },   "forprimes undef");
 ok(!eval { forprimes { 1 } 2, undef; },   "forprimes 2,undef");
@@ -62,3 +63,28 @@ ok(!eval { prime_iterator(4.5); }, "iterator 4.5");
 {my $it = prime_iterator(31398);
   is_deeply( [map { $it->() } 1..3], [31469,31477,31481], "iterator 3 primes 
starting at 31398" );
 }
+
+# For fun, nest them.
+{
+  my @t;
+  forprimes {
+    forprimes {
+      forprimes { push @t, $_ } $_,$_+10;
+    } 10*$_,10*$_+10;
+  } 10;
+  is_deeply( [@t], [qw/23 29 31 29 31 37 31 37 41 37 41 43 47 53 59 61 59 61 
67 71 73 79 73 79 83 79 83 89/], "triple nested forprimes" );
+}
+{
+  my @t;
+  my $ita = prime_iterator();
+  while ((my $a = $ita->()) <= 10) {
+    my $itb = prime_iterator(10*$a);
+    while ((my $b = $itb->()) <= 10*$a+10) {
+      my $itc = prime_iterator($b);
+      while ((my $c = $itc->()) <= $b+10) {
+        push @t, $c;
+      }
+    }
+  }
+  is_deeply( [@t], [qw/23 29 31 29 31 37 31 37 41 37 41 43 47 53 59 61 59 61 
67 71 73 79 73 79 83 79 83 89/], "triple nested iterator" );
+}
diff --git a/util.c b/util.c
index 3507eaa..f0742f4 100644
--- a/util.c
+++ b/util.c
@@ -916,6 +916,8 @@ IV _XS_mertens(UV n) {
 #endif
 }
 
+/* TODO: the two Chebyshev functions ought to use segmented sieves */
+
 double _XS_chebyshev_theta(UV n)
 {
   KAHAN_INIT(sum);

-- 
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