This is an automated email from the git hooks/post-receive script. ppm-guest pushed a commit to annotated tag v0.27 in repository libmath-prime-util-perl.
commit de8f5aa09db4c7b31667930b3ef251153a44ce8a Author: Dana Jacobsen <d...@acm.org> Date: Wed May 8 18:31:02 2013 -0700 Convenient primality proof random test --- lib/Math/Prime/Util/PP.pm | 29 ++++++++++++++----- xt/primality-proofs.pl | 72 +++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 94 insertions(+), 7 deletions(-) diff --git a/lib/Math/Prime/Util/PP.pm b/lib/Math/Prime/Util/PP.pm index a1349f7..de34252 100644 --- a/lib/Math/Prime/Util/PP.pm +++ b/lib/Math/Prime/Util/PP.pm @@ -1875,12 +1875,22 @@ sub primality_proof_bls75 { my $B = $nm1->copy; # unfactored part my @factors = (2); croak "BLS75 error: n-1 not even" unless $nm1->is_even(); - while ($B->is_even) { $B /= 2; $A *= 2; } - foreach my $f (3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79) { - last if $f*$f > $B; - if (($B % $f) == 0) { - push @factors, $f; - do { $B /= $f; $A *= $f; } while (($B % $f) == 0); + { + while ($B->is_even) { $B /= 2; $A *= 2; } + my $tlim = $_primes_small[-1]; + if ($B > $tlim*$tlim) { + my @small_primes = @_primes_small[2 .. $#_primes_small]; + foreach my $f (@small_primes) { + if (($B % $f) == 0) { + push @factors, $f; + do { $B /= $f; $A *= $f; } while (($B % $f) == 0); + last if $B <= $tlim*$tlim; + } + } + } + if ($B <= $tlim*$tlim) { + push @factors, factor($B); + $A *= $B; $B /= $B; } } my @nstack; @@ -1893,10 +1903,15 @@ sub primality_proof_bls75 { } else { push @nstack, $B; } + my $using_theorem_7 = 0; while (@nstack) { my ($s,$r) = $B->copy->bdiv($A->copy->bmul(2)); my $fpart = ($A+1) * (2*$A*$A + ($r-1) * $A + 1); last if $n < $fpart; + # Theorem 7.... + # my $tlim = $_primes_small[-1]; + # $fpart = ($A*$tlim+1) * (2*$A*$A + ($r-$tlim) * $A + 1); + # if ($n < $fpart) { $using_theorem_7 = 1; last; } my $m = pop @nstack; # Don't use bignum if it has gotten small enough. @@ -1905,8 +1920,8 @@ sub primality_proof_bls75 { my @ftry; $_holf_r = 1; foreach my $sub (@_fsublist) { - last if scalar @ftry >= 2; @ftry = $sub->($m); + last if scalar @ftry >= 2; } # If we couldn't find a factor, skip it. next unless scalar @ftry > 1; diff --git a/xt/primality-proofs.pl b/xt/primality-proofs.pl new file mode 100644 index 0000000..580be8d --- /dev/null +++ b/xt/primality-proofs.pl @@ -0,0 +1,72 @@ +#!/usr/bin/perl +use warnings; +use strict; +use Math::Prime::Util ':all'; +use Math::BigInt lib=>"GMP"; +#use Crypt::Primes 'maurer'; +use Math::Pari qw/isprime/; +use Crypt::Random 'makerandom'; +use Data::Dump 'dump'; + +$|++; +my $num = 71; +my $size = 300; +my $prime_method = 'pari'; # mpu, pari, or cpmaurer + +my @ns; +print "Generate "; +die "invalid size, must be > 4" unless $size > 4; +foreach my $i (1..$num) { + my $bits = int(rand($size-4)) + 4; + my $n; + + # How do we get random primes? + # MPU is the fastest, but it's our own code with identical primality tests. + # Pari + Crypt::Random works pretty well if you have them. + # Crypt::Primes::maurer will sometimes output composites (!!!). + if ($prime_method eq 'cpmaurer') { + $n = Crypt::Primes::maurer(Size=>$bits); + } elsif ($prime_method eq 'pari') { + do { $n = makerandom(Size=>$bits,Strength=>0); } while !isprime($n); + } elsif ($prime_method eq 'mpu') { + $n = random_nbit_prime($bits); + } else { + die "Unknown random prime generation method\n"; + } + + push @ns, Math::BigInt->new("$n"); + # print a number corresponding to hundreds of bits + print int(3.322*length("$n")/100); +} +print "\n"; + +my @certs; +print "Prove "; +foreach my $n (@ns) { + my ($isp,$cert) = is_provable_prime_with_cert($n); + die "$n is reported as $isp\n" unless $isp == 2; + push @certs, [$n, $cert]; + print proof_mark($cert); +} +print "\n"; + +print "Verify "; +foreach my $certn (@certs) { + my $v = verify_prime($certn->[1]); + print proof_mark($certn->[1]); + next if $v; + print "\n\n$certn->[0] didn't verify!\n\n"; + print dump($certn->[1]); + die; +} +print "\n"; + +sub proof_mark { + my $cert = shift; + my $type = (scalar @$cert == 1) ? "bpsw" : $cert->[1]; + if (!defined $type) { die "\nNo cert:\n\n", dump($cert); } + if ($type =~ /n-1/i) { return '1'; } + elsif ($type =~ /bpsw/i) { return '.'; } + elsif ($type =~ /ecpp|agkm/i) { return 'E'; } + return '?'; +} -- 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