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