[libmath-prime-util-perl] 03/59: Add some more implementations

Partha P. Mukherjee ppm-guest at moszumanska.debian.org
Thu May 21 18:44:51 UTC 2015


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

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

commit b728d6b7af7a5d6cfed99044d2d621a19333a49e
Author: Dana Jacobsen <dana at acm.org>
Date:   Wed Jun 27 01:22:17 2012 -0600

    Add some more implementations
---
 examples/bench-sieve-pp.pl | 215 ++++++++++++++++++++++++++++++++++++++-------
 1 file changed, 182 insertions(+), 33 deletions(-)

diff --git a/examples/bench-sieve-pp.pl b/examples/bench-sieve-pp.pl
index 6eb84d8..724351b 100644
--- a/examples/bench-sieve-pp.pl
+++ b/examples/bench-sieve-pp.pl
@@ -3,28 +3,46 @@ use strict;
 use warnings;
 
 use Benchmark qw/:all/;
+#use Devel::Size qw/total_size/;
+#use Math::Prime::Util;
+#use Math::Prime::FastSieve;
+#*mpu_erat = \&Math::Prime::Util::erat_primes;
+#*fs_erat = \&Math::Prime::FastSieve::primes;
 
 my $upper = shift || 8192;
 my $count = shift || -1;
 my $countarg;
 
+# Shows sizes for sieving to 100k, and rate/second for sieving to 16k
 my $pc_subs = {
-  "Scriptol"    => sub { scriptol($countarg) },
-  "Shootout"    => sub { shootout($countarg) },
-  "In Many"     => sub { inmany($countarg) },
-  "Rosetta 1"   => sub { rosetta1($countarg) },
-  "Rosetta 2"   => sub { rosetta2($countarg) },
-  "Rosetta 3"   => sub { rosetta3($countarg) },
-  "Rosetta 4"   => sub { rosetta4($countarg) },
-  "DJ Vec"      => sub { dj1($countarg) },
-  "DJ Array"    => sub { dj2($countarg) },
-  "DJ String1"  => sub { dj3($countarg) },
-  "DJ String2"  => sub { dj4($countarg) },
+  "Scriptol"    => sub { scriptol($countarg) },    # 3200k    290/s
+  "Shootout"    => sub { shootout($countarg) },    # 3200k    231/s
+  "In Many"     => sub { inmany($countarg) },      # 2018k    666/s
+  "Rosetta 1"   => sub { rosetta1($countarg) },    # 3449k    187/s
+  "Rosetta 2"   => sub { rosetta2($countarg) },    #   13k    109/s
+  "Rosetta 3"   => sub { rosetta3($countarg) },    # 4496k    165/s
+  "Rosetta 4"   => sub { rosetta4($countarg) },    #   25k     60/s
+  "Atkin MPTA"  => sub { atkin($countarg) },       # 3430k     90/s
+  "DO Array"    => sub {daoswald_array($countarg)},# 3200k    306/s
+  "DO Vec"      => sub {daoswald_vec($countarg)},  #   13k    112/s
+  "Merlyn"      => sub { merlyn($countarg)},       #   13k     96/s
+  "DJ Vec"      => sub { dj1($countarg) },         #    7k    245/s
+  "DJ Array"    => sub { dj2($countarg) },         # 1494k    475/s
+  "DJ String1"  => sub { dj3($countarg) },         #   50k    981/s
+  "DJ String2"  => sub { dj4($countarg) },         #   50k   1682/s
+#  "MPU Sieve"   => sub {
+#               scalar @{mpu_erat(2,$countarg)}; }, #    3k  14325/s
+#  "MPFS Sieve"  => sub {
+#               scalar @{fs_erat($countarg)};    }, #    7k  14325/s
 };
 
 my %verify = (
       10 => 4,
+      11 => 5,
      100 => 25,
+     112 => 29,
+     113 => 30,
+     114 => 30,
     1000 => 168,
    10000 => 1229,
   100000 => 9592,
@@ -48,8 +66,8 @@ cmpthese($count, $pc_subs);
 # www.scriptol.com/programming/sieve.php
 sub scriptol {
   my($max) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $max < 2;
+  return 1 if $max < 3;
 
   my @flags = (0 .. $max);
   for my $i (2 .. int(sqrt($max)) + 1)
@@ -60,6 +78,7 @@ sub scriptol {
            undef $flags[$k];
        }
   }
+  #print "scriptol size: ", total_size(\@flags), "\n" if $max > 90000;
   my $count = 0;
   for my $j (2 .. $max) {
     $count++ if defined $flags[$j];
@@ -70,8 +89,8 @@ sub scriptol {
 # http://dada.perl.it/shootout/sieve.perl.html
 sub shootout {
   my($max) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $max < 2;
+  return 1 if $max < 3;
 
   my $count = 0;
   my @flags = (0 .. $max);
@@ -82,14 +101,16 @@ sub shootout {
        }
        $count++;
   }
+  #print "shootout size: ", total_size(\@flags), "\n" if $max > 90000;
   $count;
 }
 
 # http://c2.com/cgi/wiki?SieveOfEratosthenesInManyProgrammingLanguages
 sub inmany {
   my($max) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $max < 2;
+  return 1 if $max < 3;
+  $max++;
 
   my @c;
   for(my $t=3; $t*$t<$max; $t+=2) {
@@ -97,6 +118,7 @@ sub inmany {
          for(my $s=$t*$t; $s<$max; $s+=$t*2) { $c[$s]++ }
      }
  }
+ #print "inmany size: ", total_size(\@c), "\n" if $max > 90000;
  my $count = 1;
  for(my $t=3; $t<$max; $t+=2) {
      $c[$t] || $count++;
@@ -107,8 +129,8 @@ sub inmany {
 # http://rosettacode.org/wiki/Sieve_of_Eratosthenes#Perl
 sub rosetta1 {
   my($max) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $max < 2;
+  return 1 if $max < 3;
 
   my $count = 0;  #my @primes;
   my @tested = (1);
@@ -120,14 +142,15 @@ sub rosetta1 {
       $tested[$k-1]= 1;
     }
   }
+  #print "R1 size: ", total_size(\@tested), "\n" if $max > 90000;
   $count;  #scalar @primes;
 }
 
 # http://rosettacode.org/wiki/Sieve_of_Eratosthenes#Perl
 sub rosetta2 {
   my($max) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $max < 2;
+  return 1 if $max < 3;
 
   my $count = 0; #my @primes;
   my $nonPrimes = '';
@@ -139,14 +162,15 @@ sub rosetta2 {
       $count++; #push @primes, $p;
     }
   }
+  #print "R2 size: ", total_size(\$nonPrimes), "\n" if $max > 90000;
   $count; #scalar @primes;
 }
 
 # http://rosettacode.org/wiki/Sieve_of_Eratosthenes#Perl
 sub rosetta3 {
   my($max) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $max < 2;
+  return 1 if $max < 3;
 
   my $i;
   my @s;
@@ -154,14 +178,15 @@ sub rosetta3 {
               grep { not $s[ $i  = $_ ] and do
 		 { $s[ $i += $_ ]++ while $i <= $max; 1 }
 	} 2 .. $max;
+  #print "R3 size: ", total_size(\@s), "\n" if $max > 90000;
   $count; #scalar @primes;
 }
 
 # http://rosettacode.org/wiki/Sieve_of_Eratosthenes#Perl
 sub rosetta4 {
   my($max) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $max < 2;
+  return 1 if $max < 3;
 
   my $i;
   my $s = '';
@@ -169,13 +194,133 @@ sub rosetta4 {
               grep { not vec $s, $i  = $_, 1 and do
 		 { (vec $s, $i += $_, 1) = 1 while $i <= $max; 1 }
 	} 2 .. $max;
+  #print "R4 size: ", total_size(\$s), "\n" if $max > 90000;
   $count; #scalar @primes;
 }
 
+# From Math::Primes::TiedArray
+sub atkin {
+  my($max) = @_;
+  return 0 if $max < 2;
+  return 1 if $max < 3;
+  return 2 if $max < 5;
+
+    my $sqrt     = sqrt($max);
+    my %sieve;
+    foreach my $x ( 1 .. $sqrt ) {
+
+        foreach my $y ( 1 .. $sqrt ) {
+
+            my $n = 3 * $x**2 - $y**2;
+            if (    $x > $y
+                and $n <= $max
+                and $n % 12 == 11 )
+            {
+                $sieve{$n} = not $sieve{$n};
+            }
+
+            $n = 3 * $x**2 + $y**2;
+            if ( $n <= $max and $n % 12 == 7 ) {
+                $sieve{$n} = not $sieve{$n};
+            }
+
+            $n = 4 * $x**2 + $y**2;
+            if (    $n <= $max
+                and ( $n % 12 == 1 or $n % 12 == 5 ) )
+            {
+                $sieve{$n} = not $sieve{$n};
+            }
+        }
+    }
+    # eliminate composites by sieving
+    foreach my $n ( 5 .. $sqrt ) {
+
+        next unless $sieve{$n};
+
+        my $k = int(1/$n**2) * $n**2;
+        while ( $k <= $max ) {
+            $sieve{$k} = 0;
+            $k += $n**2;
+        }
+    }
+    $sieve{2} = 1;
+    $sieve{3} = 1;
+    #print "Atkin size: ", total_size(\%sieve), "\n" if $max > 90000;
+
+    # save the found primes in our cache
+    my $count = 0;
+    foreach my $n ( 2 .. $max ) {
+        next unless $sieve{$n};
+        $count++;
+    }
+    $count;
+}
+
+# https://github.com/daoswald/Inline-C-Perl-Mongers-Talk/blob/master/primesbench.pl
+sub daoswald_array {
+  my($top) = @_;
+  return 0 if $top < 2;
+  return 1 if $top < 3;
+
+    my @primes = (1) x $top;
+    my $i_times_j;
+    for my $i ( 2 .. sqrt $top ) {
+        if ( $primes[$i] ) {
+            for ( my $j = $i; ( $i_times_j = $i * $j ) < $top; $j++ ) {
+                undef $primes[$i_times_j];
+            }
+        }
+    }
+  #print "do_array size: ", total_size(\@primes), "\n" if $top > 90000;
+  my $count = scalar grep { $primes[$_] } 2 .. $#primes;
+  $count;
+}
+
+sub daoswald_vec {
+  my($top) = @_;
+  return 0 if $top < 2;
+  return 1 if $top < 3;
+
+  my $primes = '';
+    vec( $primes, $top, 1 ) = 0;
+    my $i_times_j;
+    for my $i ( 2 .. sqrt $top ) {
+        if ( !vec( $primes, $i, 1 ) ) {
+            for ( my $j = $i; ( $i_times_j = $i * $j ) < $top; $j++ ) {
+                vec( $primes, $i_times_j, 1 ) = 1;
+            }
+        }
+    }
+  #print "do_vec size: ", total_size(\$primes), "\n" if $top > 90000;
+  my $count = scalar grep { !vec( $primes, $_, 1 ) } 2 .. $top-1 ;
+  $count;
+}
+
+# Merlyn's Unix Review Column 26, June 1999
+# http://www.stonehenge.com/merlyn/UnixReview/col26.html
+sub merlyn {
+  my($UPPER) = @_;
+  return 0 if $UPPER < 2;
+  return 1 if $UPPER < 3;
+
+  my $count = 0;
+  my $sieve = "";
+  GUESS: for (my $guess = 2; $guess <= $UPPER; $guess++) {
+    next GUESS if vec($sieve,$guess,1);
+    $count++;
+    for (my $mults = $guess * $guess; $mults <= $UPPER; $mults += $guess) {
+      vec($sieve,$mults,1) = 1;
+    }
+  }
+  #print "Merlyn size: ", total_size(\$sieve), "\n" if $UPPER > 90000;
+  $count;
+}
+
+
 sub dj1 {
   my($end) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $end < 2;
+  return 1 if $end < 3;
 
   # vector
   my $sieve = '';
@@ -188,6 +333,7 @@ sub dj1 {
     }
     do { $n += 2 } while vec($sieve, $n >> 1, 1) != 0;
   }
+  #print "DJ1 size: ", total_size(\$sieve), "\n" if $end > 90000;
   my $count = 1;
   $n = 3;
   while ($n <= $end) {
@@ -199,8 +345,8 @@ sub dj1 {
 
 sub dj2 {
   my($end) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $end < 2;
+  return 1 if $end < 3;
 
   # array
   my @sieve;
@@ -213,6 +359,7 @@ sub dj2 {
     }
     do { $n += 2 } while $sieve[$n>>1];
   }
+  #print "DJ2 size: ", total_size(\@sieve), "\n" if $end > 90000;
   my $count = 1;
   $n = 3;
   while ($n <= $end) {
@@ -226,8 +373,8 @@ sub dj2 {
 # which is just this code with a presieve added.
 sub dj3 {
   my($end) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $end < 2;
+  return 1 if $end < 3;
   $end-- if ($end & 1) == 0;
 
   # string
@@ -243,6 +390,7 @@ sub dj3 {
     }
     do { $n += 2 } while substr($sieve, $n>>1, 1);
   }
+  #print "DJ3 size: ", total_size(\$sieve), "\n" if $end > 90000;
   my $count = 1 + $sieve =~ tr/0//;
   $count;
 }
@@ -250,8 +398,8 @@ sub dj3 {
 # 2-3x faster than inmany, 6-7x faster than any of the other non-DJ methods.
 sub dj4 {
   my($end) = @_;
-  return 0 if $upper < 2;
-  return 1 if $upper < 3;
+  return 0 if $end < 2;
+  return 1 if $end < 3;
   $end-- if ($end & 1) == 0;
 
   # string with prefill
@@ -269,6 +417,7 @@ sub dj4 {
     }
     do { $n += 2 } while substr($sieve, $n>>1, 1);
   }
+  #print "DJ4 size: ", total_size(\$sieve), "\n" if $end > 90000;
   my $count = 1 + $sieve =~ tr/0//;
   $count;
 }

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