This is an automated email from the git hooks/post-receive script.
ppm-guest pushed a commit to annotated tag v0.37
in repository libmath-prime-util-perl.
commit 1f31843759f3d21a0ff9926d2bafcf1eae4dfac0
Author: Dana Jacobsen <d...@acm.org>
Date: Fri Jan 24 17:04:58 2014 -0800
Update some examples
---
examples/abundant.pl | 2 +-
examples/sophie_germain.pl | 133 ++++++++++++++++++++++++++++-----------------
examples/twin_primes.pl | 30 ++++++----
3 files changed, 103 insertions(+), 62 deletions(-)
diff --git a/examples/abundant.pl b/examples/abundant.pl
index 23b3564..8cfe2b4 100755
--- a/examples/abundant.pl
+++ b/examples/abundant.pl
@@ -5,7 +5,6 @@ use warnings;
# Find the first N abundant, deficient, or perfect numbers.
use Math::Prime::Util qw/divisor_sum next_prime is_prime/;
-use Math::BigInt try => "GMP,Pari";
my $count = shift || 20;
my $type = lc(shift || 'abundant');
@@ -26,6 +25,7 @@ if ($type eq 'abundant') {
# We just look for 2^(p-1)*(2^p-1) where 2^p-1 is prime.
# Basically we're just finding Mersenne primes.
# It's possible there are odd perfect numbers larger than 10^1500.
+ do { require Math::BigInt; Math::BigInt->import(try=>"GMP,Pari"); };
while ($count-- > 0) {
while (1) {
$p = next_prime($p);
diff --git a/examples/sophie_germain.pl b/examples/sophie_germain.pl
index 0fec4b1..971ad3b 100755
--- a/examples/sophie_germain.pl
+++ b/examples/sophie_germain.pl
@@ -6,59 +6,92 @@ use Math::Prime::Util qw/prime_iterator is_prime
next_prime nth_prime_upper prime_precalc forprimes/;
my $count = shift || 20;
+my $method = shift || 'forprimes';
+my $precalc = 0; # If set, precalc all the values we'll call is_prime on
# Find Sophie Germain primes (numbers where p and 2p+1 are both prime).
-# In this method, we add a filter in front of our iterator, to create a
-# Sophie-Germain-prime iterator. This isn't the fastest way, but it's still
-# 20x faster than Math::NumSeq::SophieGermainPrimes at 300k. If we add the
-# two-line precalc shown below, we can get another 4x or more.
-#
-# Example:
-# time perl examples/sophie_germain.pl 300000 | md5sum
-# d380d31256cc9bc54eb5f236b3edc16d -
-# 9.673s
-#
-# time perl -MMath::NumSeq::SophieGermainPrimes -E 'my $seq =
Math::NumSeq::SophieGermainPrimes->new; do { say 0+($seq->next)[1] } for
1..300000' | md5sum
-# d380d31256cc9bc54eb5f236b3edc16d -
-# 4m11.5s
+# Four methods are shown: forprimes, iter, iter2, and MNS.
+
+# Times for 300k:
#
-# With method 2:
-# time perl examples/sophie_germain.pl 300000 | md5sum
-# d380d31256cc9bc54eb5f236b3edc16d -
-# 1.828s
-
-
-sub get_sophie_germain_iterator {
- my $p = shift || 2;
- my $it = prime_iterator($p);
- return sub {
- do { $p = $it->() } while !is_prime(2*$p+1);
- $p;
- };
+# 300k 1M
+# precalc:
+# forprimes 1.3s 9.0MB 7.1s 21.6MB
+# iter 2.8s 8.7MB 12.6s 21.4MB
+# iter2 1.9s 8.7MB 9.4s 21.4MB
+# no precalc:
+# forprimes 1.5s 4.5MB 5.6s 4.5MB
+# iter 9.5s 4.3MB 37.5s 4.3MB
+# iter2 8.5s 4.3MB 33.9s 4.3MB
+# MNS 254.3s 11.3MB >1500s >15 MB
+
+if ($precalc) {
+ prime_precalc(2 * sg_upper_bound($count));
}
-my $sgit = get_sophie_germain_iterator();
-print $sgit->(), "\n" for 1 .. $count;
-# Method 2. At 300k this is 70x faster than Math::NumSeq::SophieGermainPrimes.
-#
-#my $estimate = 100 + int( nth_prime_upper($count) * 1.6 * log($count) );
-#prime_precalc(2 * $estimate);
-#
-#my $prime = 2;
-#for (1..$count) {
-# $prime = next_prime($prime) while (!is_prime(2*$prime+1));
-# print "$prime\n";
-# $prime = next_prime($prime);
-#}
-
-# Alternate method, 10-20% faster, would benefit from a tighter estimate.
-#
-# my $numfound = 0;
-# forprimes {
-# if ($numfound < $count && is_prime(2*$_+1)) {
-# print "$_\n";
-# $numfound++;
-# }
-# } $estimate;
-# die "Estimate too low" unless $numfound >= $count;
+
+if ($method eq 'forprimes') {
+
+ my $estimate = sg_upper_bound($count);
+ my $numfound = 0;
+ forprimes {
+ if ($numfound < $count && is_prime(2*$_+1)) {
+ print "$_\n";
+ $numfound++;
+ }
+ } $estimate;
+ die "Estimate too low" unless $numfound >= $count;
+
+} elsif ($method eq 'iter') {
+
+ # Wrap the standard iterator
+ sub get_sophie_germain_iterator {
+ my $p = shift || 2;
+ my $it = prime_iterator($p);
+ return sub {
+ do { $p = $it->() } while !is_prime(2*$p+1);
+ $p;
+ };
+ }
+ my $sgit = get_sophie_germain_iterator();
+ print $sgit->(), "\n" for 1 .. $count;
+
+} elsif ($method eq 'iter2') {
+
+ # Iterate directly using next_prime
+ my $prime = 2;
+ for (1 .. $count) {
+ $prime = next_prime($prime) while !is_prime(2*$prime+1);
+ print "$prime\n";
+ $prime = next_prime($prime);
+ }
+
+} elsif ($method eq 'MNS') {
+
+ # Use Math::NumSeq
+ require Math::NumSeq::SophieGermainPrimes;
+ my $seq = Math::NumSeq::SophieGermainPrimes->new;
+ for (1 .. $count) {
+ print 0+($seq->next)[1];
+ }
+
+}
+
+# Used for precalc and the forprimes example
+sub sg_upper_bound {
+ my $count = shift;
+ my $nth = nth_prime_upper($count);
+ # For lack of a better formula, do this step-wise estimate.
+ my $estimate = ($count < 5000) ? 150 + int( $nth * log($nth) * 1.2 )
+ : ($count < 19000) ? int( $nth * log($nth) * 1.135 )
+ : ($count < 45000) ? int( $nth * log($nth) * 1.10 )
+ : ($count < 100000) ? int( $nth * log($nth) * 1.08 )
+ : ($count < 165000) ? int( $nth * log($nth) * 1.06 )
+ : ($count < 360000) ? int( $nth * log($nth) * 1.05 )
+ : ($count < 750000) ? int( $nth * log($nth) * 1.04 )
+ : ($count <1700000) ? int( $nth * log($nth) * 1.03 )
+ : int( $nth * log($nth) * 1.02 );
+
+ return $estimate;
+}
diff --git a/examples/twin_primes.pl b/examples/twin_primes.pl
index 362a737..514cc70 100755
--- a/examples/twin_primes.pl
+++ b/examples/twin_primes.pl
@@ -2,8 +2,7 @@
use strict;
use warnings;
-use Math::Prime::Util qw/-nobigint
- prime_iterator prime_iterator_object
+use Math::Prime::Util qw/prime_iterator prime_iterator_object
next_prime is_prime
nth_prime_upper prime_precalc/;
@@ -12,16 +11,25 @@ my $count = shift || 20;
# Find twin primes (numbers where p and p+2 are prime)
# Time for the first 300k:
+#
+# Not iterators:
+# 0.6s forprimes { say $l if $l+2==$_; $l=$_; } 64764841
# 1.0s bin/primes.pl --twin 2 64764839
-# 1.4s get_twin_prime_iterator2
-# 2.3s get_twin_prime_iterator1
-# 4.1s get_twin_prime_iterator3
-# 6.1s get_twin_prime_iterator3 (object iterator)
-# 7.6s get_twin_prime_iterator2 without precalc
-# 8.4s get_twin_prime_iterator1 without precalc
-# 10.9s get_twin_prime_iterator3 without precalc
-# 13.1s get_twin_prime_iterator3 without precalc (object iterator)
-# 219.8s Math::NumSeq::TwinPrimes (Perl 5.19.4 with v66)
+#
+# Iterators with precalc:
+# 1.6s get_twin_prime_iterator2
+# 2.4s get_twin_prime_iterator1
+# 4.2s get_twin_prime_iterator3
+# 4.5s get_twin_prime_iterator4 (object iterator)
+#
+# Iterators without precalc:
+# 7.7s get_twin_prime_iterator2
+# 8.5s get_twin_prime_iterator1
+# 10.8s get_twin_prime_iterator3
+# 16.7s get_twin_prime_iterator4 (object iterator)
+#
+# Alternatives:
+# 251.9s Math::NumSeq::TwinPrimes (Perl 5.19.7, Math::NumSeq 67)
# This speeds things up, but isn't necessary.
my $estimate = 5000 + int( nth_prime_upper($count) * 1.4 * log($count) );
--
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
