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