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 9a3fc47eb791c223c03737ac44b61388b719ad4c
Author: Dana Jacobsen <d...@acm.org>
Date: Fri Jan 24 18:35:19 2014 -0800
Update examples and examples/README
---
examples/README | 119 +++++++++++-------------------------------
examples/find_mr_bases.pl | 46 +++++++---------
examples/parallel_fibprime.pl | 18 +++++--
examples/verify-cert.pl | 35 +------------
4 files changed, 66 insertions(+), 152 deletions(-)
diff --git a/examples/README b/examples/README
index ac51401..1a40ac1 100644
--- a/examples/README
+++ b/examples/README
@@ -1,108 +1,51 @@
-There are two main types of scripts here: benchmarks and correctness tests.
+abundant.pl
-The test-* scripts are generally trying to test one part of the module
-against another part of the module, an external module, or an external program.
-These usually consist of a combination of fixed tests and a long sequence of
-testing with random numbers, trying to find things the standard testing might
-have missed.
+ Prints the first N abundant (or deficient, or perfect) numbers. E.g:
+ perl abundant.pl 100 abundant
+ perl abundant.pl 100 deficient
+ perl abundant.pl 15 perfect
-test-factor-yafu.pl
+sophie_germain.pl
- Tests factorization compared with YAFU (v1.31.1). No arguments.
+ Prints the first N Sophie-Germain primes. E.g.:
-test-factor-mpxs.pl
+ perl sophia_germain.pl 100000
- Tests factorization compared with Math::Factor::XS (v0.26).
- One argument gives the number of random tests to perform.
+twin_primes.pl
-test-nextprime-yafu.pl
+ Prints the first N twin-primes (first value of the pair). E.g.:
- Tests next_prime() compared with YAFU (v1.31.1). No arguments.
+ perl twin_primes.pl 100000
-test-primes-yafu.pl
- Tests primes($a,$b+$interval) compared with YAFU (v1.31.1). No arguments.
- The interval is currently 8000.
+find_mr_bases.pl
+
+ An example using threads to do a parallel search for good deterministic
+ bases for a Miller-Rabin test. This is definitely not the fastest way
+ to find these, but it's a decent example of quickly trying out an idea.
+ Be sure to set $nthreads to the right value for your machine. It should
+ fully load your CPUs.
-test-holf.pl
- Tests the holf_factor() function vs. the factor() function. Given enough
- rounds, HOLF (like Fermat) should be able to factor a number. We keep
- calling it on each non-prime return value until it's done.
+parallel_fibprime.pl
-test-nthapprox.pl
+ Find Fibonacci primes, in parallel. You will want Math::Prime::Util::GMP
+ installed, as these are many-thousand-digit numbers.
- Tests the nth_prime approximation and upper/lower bounds vs. known values
- for the nth prime on large values.
+verify-cert.pl
-test-pcapprox.pl
+ Takes an MPU or Primo primality certificate and verifies it. This is
+ obsolete, as Math::Prime::Util::GMP now includes C code for this.
- Tests the prime_count approximation and upper/lower bounds vs. known values
- for Pi(x) on large values. Also examines the Schoenfeld and Stoll
- inequalities.
+verify-gmp-ecpp-cert.pl
-bench-factor-extra.pl
+ Parses the verbose output of GMP-ECPP to construct a certificate, then
+ runs it through the verification process.
- Benchmarks the various factoring options (prho, pbrent, pminus1, fermat,
holf,
- squfof, trial) on random n-digit numbers. Also gives the percent of the time
- solutions were found with the given number of rounds (256k for SQUFOF, 2000
- for HOLF, and 400 for probabilistic methods).
+verify-sage-ecpp-cert.pl
-bench-factor.pl
-
- Benchmarks factoring random and semiprime n-digit numbers using the factor
- method, and compares vs. Math::Factor::XS (v0.26). The latter uses a trial
- division algorithm. MPU 0.05 and later use a threshold of 10M (8 digits)
- to switch between trial division and methods like SQUFOF, Pollard's Rho, and
- HOLF.
-
-bench-factor-semiprime.pl
-
- Benchmarks factoring semiprimes, and compares with Math::Prime::XS and
- Math::Pari. Takes two optional arguments: the number of digits (default 15)
- and the benchmark count (default -2, meaning 2 seconds).
-
-bench-is-prime.pl
-
- Benchmarks is_prime on random n-digit numbers, n from 5 to 10/20. Compares
- MPU's is_prime and is_prob_prime vs. Math::Prime::XS and optionally
- Math::Primality, Math::Pari, and Math::Prime::FastSieve. The first two of
- the optional modules use methods more appropriate for big numbers, so are up
- to an order of magnitude slower for 64-bit numbers. The last module (MPFS)
- is extremely fast, but requires presieving to at least the number to be
- tested, which is great for small numbers, but not for large.
- Also, no additional precalc is done for MPU. If you really want blazing
- fast is_prime and don't care about the memory and time to sieve, run
- prime_precalc to the limit of your numbers and is_prime will turn into a
- sub-100 microsecond bit array lookup for any number in the range.
- Takes one optional argument of the benchmark count (default -5).
-
-bench-miller-rabin.pl
-
- Benchmarks the strong miller_rabin test using 7 bases at various digit
counts.
- Takes one optional argument of the benchmark count (default -5).
-
-bench-nthprime.pl
-
- Benchmarks the speed of nth_prime with various digit sizes.
-
-bench-pcapprox.pl
-
- Benchmarks the speed of prime_count related functions for random n-digit
- numbers (n = 5 to 10/20). This includes lower bound, lower/upper average,
- the prime_count_approx function, li(x), and R(x).
- Takes one optional argument of the benchmark count (default -5).
-
-bench-primecount.pl
-
- Benchmarks the speed of prime_count on random n-digit numbers (n = 2 .. 8).
- Takes one optional argument of the benchmark count (default -5).
- I'll note you can easily see the transition from where we're just counting
- existing value to where we have to sieve. Adding a prime_precalc(10**9)
- line will speed up the 5-,6-,7-, and 8-digit prime_counts greatly.
-
-bench-random-prime.pl
-
- Benchmarks the speed of random_ndigit_prime for various digits.
+ Verifies the output of SAGE's ECPP. The SAGE module looks like it died
+ in development and never got into SAGE. NZMath's ECPP doesn't seem to
+ output a certificate, which makes it much less useful.
diff --git a/examples/find_mr_bases.pl b/examples/find_mr_bases.pl
index 8f78173..47ba7d1 100755
--- a/examples/find_mr_bases.pl
+++ b/examples/find_mr_bases.pl
@@ -3,15 +3,13 @@ use warnings;
use strict;
use threads;
use threads::shared;
-use Math::Prime::Util qw/is_prime is_strong_pseudoprime/;
-my $nthreads = 12;
+use Math::Prime::Util qw/is_prime is_strong_pseudoprime forcomposites/;
+my $nthreads = 4;
# Single base.
my @composites;
-for my $n (3 .. 1000000) {
- push @composites, $n if $n % 2 && !is_prime($n);
-}
+forcomposites { push @composites, $_ if $_ % 2; } 1_000_000;
# Serial:
# my $base = 2;
@@ -31,52 +29,45 @@ for my $n (3 .. 1000000) {
# Parallel:
my $maxn :shared;
-my $start = int(2**59+2**41); # People have mined below 2^55
+my $start = int(2**60+2**41); # People have mined below 2^55
$maxn = 2047;
-my $nextn = 2049;
my @threads;
-push @threads, threads->create('search_bases', $start, $_) for
(0..$nthreads-1);
+push @threads, threads->create('search_bases', $start, $_) for 1..$nthreads;
# We should sit here doing cond_waits on a results array.
$_->join() for (@threads);
sub search_bases {
my($start, $t) = @_;
- my $base = $start + $t;
- while (1) {
- do { $base += $t; next; } if is_strong_pseudoprime($nextn, $base);
+ for (my $base = $start + $t - 1; 1; $base += $t) {
+ next if is_strong_pseudoprime(4, $base) || is_strong_pseudoprime(6, $base);
for my $n (@composites) {
if (is_strong_pseudoprime($n,$base)) {
if ($n > $maxn) {
lock($maxn);
- print "base $base good up to $n\n";
+ print "base $base good up to $n\n" if $n > $maxn;
$maxn = $n;
- $nextn = $n+2; $nextn++ while is_prime($nextn);
}
last;
}
}
- $base += $t;
}
}
__END__
base 2 good up to 2047
-base 1320 good up to 4097
-base 4712 good up to 4711
-base 5628 good up to 5627
-base 7252 good up to 7251
-base 7852 good up to 7851
-base 14787 good up to 9409
-base 17340 good up to 10261
-base 61380 good up to 11359
-base 78750 good up to 13747
+base 3273 good up to 2209
+base 4414 good up to 2443
+base 5222 good up to 2611
+base 8286 good up to 4033
+base 10822 good up to 5411
+base 13011 good up to 6505
+base 67910 good up to 9073
+base 82967 good up to 10371
base 254923 good up to 18299
-base 486605 good up to 25761
-base 1804842 good up to 32761
+base 2974927 good up to 18721
base 4095086 good up to 38323
-base 12772344 good up to 40501
-base 42162995 good up to 97921
+base 70903283 good up to 38503
(best results known, not found with this program)
2011-02-12 base 814494960528 good up to 132239
@@ -84,3 +75,4 @@ base 42162995 good up to 97921
2012-10-15 base 1769236083487960 good up to 192001
2012-10-17 base 1948244569546278 good up to 212321
2013-01-14 base 34933608779780163 good up to 218245
+2013-03-03 base 9345883071009581737 good up to 341531
diff --git a/examples/parallel_fibprime.pl b/examples/parallel_fibprime.pl
index 138d3e5..e9a6361 100755
--- a/examples/parallel_fibprime.pl
+++ b/examples/parallel_fibprime.pl
@@ -3,7 +3,19 @@ use strict;
use warnings;
use threads;
use threads::shared;
-use Math::BigInt lib => 'GMP';
+
+# Overkill, but let's try to select a good bigint module.
+my $bigint_class;
+if (eval { require Math::GMPz; 1; }) {
+ $bigint_class = "Math::GMPz";
+} elsif (eval { require Math::GMP; 1; }) {
+ $bigint_class = "Math::GMP";
+} else {
+ require Math::BigInt;
+ Math::BigInt->import(try=>"GMP,Pari");
+ $bigint_class = "Math::BigInt";
+}
+
use Math::Prime::Util ':all';
use Time::HiRes qw(gettimeofday tv_interval);
$| = 1;
@@ -51,7 +63,7 @@ my @found :shared; # push the primes found here
my @karray : shared; # array of min k for each thread
my @threads;
-push @threads, threads->create('fibprime', $_) for (1..$nthreads);
+push @threads, threads->create('fibprime', $_) for 1 .. $nthreads;
# Let the threads work for a little before starting the display loop
sleep 2;
@@ -80,7 +92,7 @@ $_->join() for (@threads);
sub fib_n {
my ($n, $fibstate) = @_;
- @$fibstate = (1, Math::BigInt->new(0), Math::BigInt->new(1))
+ @$fibstate = (1, $bigint_class->new(0), $bigint_class->new(1))
unless defined $fibstate->[0];
my ($curn, $a, $b) = @$fibstate;
die "fib_n only increases" if $n < $curn;
diff --git a/examples/verify-cert.pl b/examples/verify-cert.pl
index 98ead86..70d8b2d 100755
--- a/examples/verify-cert.pl
+++ b/examples/verify-cert.pl
@@ -367,7 +367,7 @@ sub verify_bls15 {
fail "BLS15: $n failed 2Q-1 > sqrt(N)" unless 2*$q-1 > $n->copy->bsqrt();
my $D = $lp*$lp - 4*$lq;
fail "BLS15: $n failed D != 0" unless $D != 0;
- fail "BLS15: $n failed jacobi(D,N) = -1" unless _jacobi($D,$n) == -1;
+ fail "BLS15: $n failed jacobi(D,N) = -1" unless kronecker($D,$n) == -1;
fail "BLS15: $n failed V_{m/2} mod N != 0"
unless (lucas_sequence($n, $lp, $lq, $m/2))[1] != 0;
fail "BLS15: $n failed V_{(N+1)/2} mod N == 0"
@@ -563,36 +563,3 @@ sub _is_perfect_square {
}
0;
}
-
-# Calculate Jacobi symbol (M|N)
-sub _jacobi {
- my($n, $m) = @_;
- return 0 if $m <= 0 || ($m % 2) == 0;
- my $j = 1;
- if ($n < 0) {
- $n = -$n;
- $j = -$j if ($m % 4) == 3;
- }
- # Split loop so we can reduce n/m to non-bigints after first iteration.
- if ($n != 0) {
- while (($n % 2) == 0) {
- $n >>= 1;
- $j = -$j if ($m % 8) == 3 || ($m % 8) == 5;
- }
- ($n, $m) = ($m, $n);
- $j = -$j if ($n % 4) == 3 && ($m % 4) == 3;
- $n = $n % $m;
- $n = int($n->bstr) if ref($n) eq 'Math::BigInt' && $n <= ''.~0;
- $m = int($m->bstr) if ref($m) eq 'Math::BigInt' && $m <= ''.~0;
- }
- while ($n != 0) {
- while (($n % 2) == 0) {
- $n >>= 1;
- $j = -$j if ($m % 8) == 3 || ($m % 8) == 5;
- }
- ($n, $m) = ($m, $n);
- $j = -$j if ($n % 4) == 3 && ($m % 4) == 3;
- $n = $n % $m;
- }
- return ($m == 1) ? $j : 0;
-}
--
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
