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 17559c6acaf9d010e87060b45374aa7a527eda95
Author: Dana Jacobsen <d...@acm.org>
Date: Fri May 10 16:35:25 2013 -0700
Allow BLS75 theorem 7 n-1 proofs
---
examples/verify-gmp-ecpp-cert.pl | 5 ++-
lib/Math/Prime/Util.pm | 65 ++++++++++++++++++++++-----
lib/Math/Prime/Util/PP.pm | 94 +++++++++++++++++++++++++---------------
xt/primality-proofs.pl | 11 +++--
4 files changed, 123 insertions(+), 52 deletions(-)
diff --git a/examples/verify-gmp-ecpp-cert.pl b/examples/verify-gmp-ecpp-cert.pl
index de3e616..7c629ce 100755
--- a/examples/verify-gmp-ecpp-cert.pl
+++ b/examples/verify-gmp-ecpp-cert.pl
@@ -4,6 +4,9 @@ use strict;
use Math::BigInt try=>"GMP,Pari";
use Math::Prime::Util qw/:all/;
use Data::Dump qw/dump/;
+my $bifilter = sub { my($ctx, $n) = @_;
+ return {dump=>"$n"} if ref($n) eq "Math::BigInt";
+ undef; };
# Takes the output of GMP-ECPP, creates a certificate in the format used
# by MPU, and runs it through the verifier.
@@ -57,5 +60,5 @@ while (<>) {
}
}
-print dump(\@cert), "\n";
+print dump_filtered(\@cert, $bifilter), "\n";
print verify_prime(@cert) ? "SUCCESS\n" : "FAILURE\n";
diff --git a/lib/Math/Prime/Util.pm b/lib/Math/Prime/Util.pm
index 0a7f465..bd59e50 100644
--- a/lib/Math/Prime/Util.pm
+++ b/lib/Math/Prime/Util.pm
@@ -234,9 +234,8 @@ sub _validate_positive_integer {
if (ref($_[0])) {
$_[0] = Math::BigInt->new("$_[0]") unless ref($_[0]) eq 'Math::BigInt';
- $_bigint_small = Math::BigInt->new("$_Config{'maxparam'}")
- unless defined $_bigint_small;
- if ($_[0]->bacmp($_bigint_small) <= 0) {
+ # Stupid workaround for Math::BigInt::GMP RT # 71548
+ if ($_[0]->bacmp(''.~0) <= 0) {
$_[0] = int($_[0]->bstr);
} else {
$_[0]->upgrade(undef) if $_[0]->upgrade(); # Stop BigFloat upgrade
@@ -251,8 +250,8 @@ sub _validate_positive_integer {
}
}
# One of these will be true:
- # 1) $n <= max and $n is not a bigint
- # 2) $n > max and $n is a bigint
+ # 1) $n <= ~0 and $n is not a bigint
+ # 2) $n > ~0 and $n is a bigint
1;
}
@@ -1802,6 +1801,13 @@ sub verify_prime {
if ($method eq 'n-1') {
# BLS75 or generalized Pocklington
#
http://www.ams.org/journals/mcom/1975-29-130/S0025-5718-1975-0384673-1/S0025-5718-1975-0384673-1.pdf
+ # Pull of optional theorem 7 data
+ my ($theorem, $t7_B1, $t7_B, $t7_a) = (5, undef, undef, undef);
+ if (scalar @pdata == 3 && ref($pdata[0]) eq 'ARRAY' && $pdata[0]->[0] =~
/^(B|T7|Theorem\s*7)$/i) {
+ $theorem = 7;
+ (undef, $t7_B1, $t7_B, $t7_a) = @{shift @pdata};
+ $t7_B = Math::BigInt->new("$t7_B");
+ }
if (scalar @pdata != 2 || (ref($pdata[0]) ne 'ARRAY') || (ref($pdata[1])
ne 'ARRAY')) {
warn "verify_prime: incorrect n-1 format, must have factors and a
values\n";
return 0;
@@ -1858,10 +1864,42 @@ sub verify_prime {
print "primality fail: A and B not coprime\n" if $verbose;
return 0;
}
- # Theorem 5, m = 1, page 624
- {
+ if ($theorem == 7) {
+ if ($B != $t7_B) {
+ print "primality fail: T7 unfactored != unfactored\n" if $verbose;
+ return 0;
+ }
+ if ($t7_B1 < 1) {
+ print "primality fail: T7 B value < 1\n" if $verbose;
+ return 0;
+ }
+ # 1. Check $B has no factors smaller than $t7_B1
+ my $no_small_factors = 0;
+ if ($_HAVE_GMP) {
+ my @trial = Math::Prime::Util::GMP::trial_factor($B, $t7_B1);
+ $no_small_factors = (scalar @trial == 1);
+ } elsif ($B <= ''.~0) {
+ my @trial = Math::Prime::Util::PP::trial_factor($B, $t7_B1);
+ $no_small_factors = (scalar @trial == 1);
+ } else {
+ # This is slow when B1 > 1M, but with a bigint B it's faster than
+ # doing trial divisions (but much slower with native B).
+ $no_small_factors =
+ (Math::BigInt::bgcd($B, primorial(Math::BigInt->new($t7_B1))) ==
1);
+ }
+ if (!$no_small_factors) {
+ print "primality fail: T7 B has a factor smaller than B1\n" if
$verbose;
+ return 0;
+ }
+ # 2. Add $B and $t7_a to lists so they get checked as part of (I).
+ push @prime_factors, $B;
+ push @pfas, $t7_a;
+ }
+ { # Theorem 5, m = 1, page 624; Theorem 7, page 626
my ($s,$r) = $B->copy->bdiv($A->copy->bmul(2));
- my $fpart = ($A+1) * (2*$A*$A + ($r-1) * $A + 1);
+ my $fpart = ($theorem == 7)
+ ? ($A*$t7_B1+1) * (2*$A*$A + ($r-$t7_B1) * $A + 1)
+ : ($A+1) * (2*$A*$A + ($r-1) * $A + 1);
if ($n >= $fpart) {
print "primality fail: not enough factors\n" if $verbose;
return 0;
@@ -1883,7 +1921,7 @@ sub verify_prime {
return 0;
}
}
- print "primality success: $n by BLS75 theorem 5\n" if $verbose > 1;
+ print "primality success: $n by BLS75 theorem $theorem\n" if $verbose > 1;
return 1;
}
@@ -2873,9 +2911,9 @@ A certificate is an array holding an C<n-cert>. An
C<n-cert> is one of:
4 a^(n-1) = 1 mod n
5 a^((n-1)/f) != 1 mod n for each factor
- n,"n-1",[n-cert, ...],[a,...]
+ n,"n-1",[optional B-block],[n-cert, ...],[a,...]
An n-1 certificate suitable for the generalized Pocklington or the
- BLS75 (Brillhart-Lehmer-Selfridge 1975, theorem 5) test. The
+ BLS75 (Brillhart-Lehmer-Selfridge 1975, theorem 5) n-1 test. The
proof is performed using BLS75 theorem 5 which requires n-1 to be
factored up to (n/2)^1/3. If n-1 is factored to more than
sqrt(n), then the conditions are identical to the generalized
@@ -2892,6 +2930,11 @@ A certificate is an array holding an C<n-cert>. An
C<n-cert> is one of:
5 for each pair (f,a) representing a factor n-cert and its 'a':
- a^(n-1) = 1 mod n
- gcd( a^((n-1)/f)-1, n ) = 1
+ If the optional B block is present, then theorem 7 will be used.
+ The B-block consists of 4 items: "B" as an identifier, the
+ factoring limit B indicating that the unfactored portion has no
+ factors smaller than B, the unfactored amount F, and an 'a' value
+ to be tested with F as in step 5.
n,"AGKM",[ec-block],[ec-block],...
An Elliptic Curve certificate. We are given n, the method "AGKM"
diff --git a/lib/Math/Prime/Util/PP.pm b/lib/Math/Prime/Util/PP.pm
index de34252..5ed8e7b 100644
--- a/lib/Math/Prime/Util/PP.pm
+++ b/lib/Math/Prime/Util/PP.pm
@@ -69,19 +69,25 @@ sub _validate_positive_integer {
if ref($n) ne 'Math::BigInt' && $n =~ tr/0123456789//c;
croak "Parameter '$n' must be >= $min" if defined $min && $n < $min;
croak "Parameter '$n' must be <= $max" if defined $max && $n > $max;
- if ($n <= ~0) {
- $_[0] = $_[0]->as_number() if ref($_[0]) eq 'Math::BigFloat';
- $_[0] = int($_[0]->bstr) if ref($_[0]) eq 'Math::BigInt';
- } elsif (ref($n) ne 'Math::BigInt') {
- croak "Parameter '$n' outside of integer range" if !defined
$bigint::VERSION;
- $_[0] = Math::BigInt->new("$n"); # Make $n a proper bigint object
- $_[0]->upgrade(undef) if $_[0]->upgrade(); # Stop BigFloat upgrade
+ if (ref($_[0])) {
+ $_[0] = Math::BigInt->new("$_[0]") unless ref($_[0]) eq 'Math::BigInt';
+ # Stupid workaround for Math::BigInt::GMP RT # 71548
+ if ($_[0]->bacmp(''.~0) <= 0) {
+ $_[0] = int($_[0]->bstr);
+ croak "Didn't convert!" if ref($_[0]);
+ } else {
+ $_[0]->upgrade(undef) if $_[0]->upgrade(); # Stop BigFloat upgrade
+ }
} else {
- $_[0]->upgrade(undef) if $_[0]->upgrade(); # Stop BigFloat upgrade
+ if ($n > ~0) {
+ croak "Parameter '$n' outside of integer range" if !defined
$bigint::VERSION;
+ $_[0] = Math::BigInt->new("$n"); # Make $n a proper bigint object
+ $_[0]->upgrade(undef) if $_[0]->upgrade(); # Stop BigFloat upgrade
+ }
}
# One of these will be true:
- # 1) $n <= max and $n is not a bigint
- # 2) $n > max and $n is a bigint
+ # 1) $n <= ~0 and $n is not a bigint
+ # 2) $n > ~0 and $n is a bigint
1;
}
@@ -1115,21 +1121,32 @@ sub _basic_factor {
}
sub trial_factor {
- my($n) = @_;
+ my($n, $maxlim) = @_;
_validate_positive_integer($n);
+ $maxlim = $n unless defined $maxlim && _validate_positive_integer($maxlim);
- my @factors = _basic_factor($n);
+ # Don't use _basic factor here -- they want a trial forced.
+ #my @factors = _basic_factor($n);
+ return ($n) if $n < 4;
+ my @factors;
+ while ( !($n % 2) ) { push @factors, 2; $n = int($n / 2); }
+ while ( !($n % 3) ) { push @factors, 3; $n = int($n / 3); }
+ while ( !($n % 5) ) { push @factors, 5; $n = int($n / 5); }
+ $n = int($n->bstr) if ref($n) eq 'Math::BigInt' && $n <= ''.~0;
return @factors if $n < 4;
my $limit = int( sqrt($n) + 0.001);
+ $limit = $maxlim if $limit > $maxlim;
my $f = 7;
SEARCH: while ($f <= $limit) {
foreach my $finc (4, 2, 4, 2, 4, 6, 2, 6) {
if ( (($n % $f) == 0) && ($f <= $limit) ) {
do { push @factors, $f; $n = int($n/$f); } while (($n % $f) == 0);
- $n = int($n->bstr) if ref($n) eq 'Math::BigInt' && $n <= ~0;
- last SEARCH if $n == 1 || Math::Prime::Util::is_prob_prime($n);
+ $n = int($n->bstr) if ref($n) eq 'Math::BigInt' && $n <= ''.~0;
+ #last SEARCH if $n == 1 || Math::Prime::Util::is_prob_prime($n);
+ last SEARCH if $n == 1;
$limit = int( sqrt($n) + 0.001);
+ $limit = $maxlim if $limit > $maxlim;
}
$f += $finc;
}
@@ -1875,22 +1892,15 @@ sub primality_proof_bls75 {
my $B = $nm1->copy; # unfactored part
my @factors = (2);
croak "BLS75 error: n-1 not even" unless $nm1->is_even();
+ my $trial_B = 20000;
{
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 @tf = trial_factor($B, $trial_B);
+ pop @tf if $tf[-1] > $trial_B;
+ foreach my $f (@tf) {
+ next if $f == $factors[-1];
+ push @factors, $f;
+ do { $B /= $f; $A *= $f; } while (($B % $f) == 0);
}
}
my @nstack;
@@ -1903,15 +1913,14 @@ sub primality_proof_bls75 {
} else {
push @nstack, $B;
}
- my $using_theorem_7 = 0;
+ my $theorem = 5;
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; }
+ $fpart = ($A*$trial_B+1) * (2*$A*$A + ($r-$trial_B) * $A + 1);
+ if ($n < $fpart) { $theorem = 7; last; }
my $m = pop @nstack;
# Don't use bignum if it has gotten small enough.
@@ -1948,7 +1957,9 @@ sub primality_proof_bls75 {
}
# Did we factor enough?
my ($s,$r) = $B->copy->bdiv($A->copy->bmul(2));
- my $fpart = ($A+1) * (2*$A*$A + ($r-1) * $A + 1);
+ my $fpart = ($theorem == 5)
+ ? ($A+1) * (2*$A*$A + ($r-1) * $A + 1)
+ : ($A*$trial_B+1) * (2*$A*$A + ($r-$trial_B) * $A + 1);
return (1,[]) if $n >= $fpart;
# Check we didn't mess up
croak "BLS75 error: $A * $B != $nm1" unless $A*$B == $nm1;
@@ -1964,7 +1975,7 @@ sub primality_proof_bls75 {
foreach my $f (@factors) {
my $success = 0;
foreach my $a (2 .. 10000) {
- my $ap = Math::BigInt->new("$a");
+ my $ap = Math::BigInt->new($a);
next unless $ap->copy->bmodpow($nm1, $n) == 1;
next unless Math::BigInt::bgcd($ap->copy->bmodpow($nm1/$f, $n)->bsub(1),
$n) == 1;
push @as, $a;
@@ -1979,8 +1990,19 @@ sub primality_proof_bls75 {
}
push @fac_proofs, (scalar @$fproof == 1) ? $fproof->[0] : $fproof;
}
- my @proof = ($n, "n-1", [@fac_proofs], [@as]);
- return (2, [@proof]);
+ if ($theorem == 5) {
+ return (2, [$n, "n-1", [@fac_proofs], [@as]]);
+ } else {
+ my $t7a = 0;
+ my $f = $B;
+ foreach my $a (2 .. 10000) {
+ my $ap = Math::BigInt->new($a);
+ next unless $ap->copy->bmodpow($nm1, $n) == 1;
+ next unless Math::BigInt::bgcd($ap->copy->bmodpow($nm1/$f, $n)->bsub(1),
$n) == 1;
+ return (2, [$n, "n-1", ["B", $trial_B, $B, $a], [@fac_proofs], [@as]]);
+ }
+ }
+ return @composite;
}
diff --git a/xt/primality-proofs.pl b/xt/primality-proofs.pl
index 580be8d..4c0cf7b 100644
--- a/xt/primality-proofs.pl
+++ b/xt/primality-proofs.pl
@@ -6,7 +6,10 @@ use Math::BigInt lib=>"GMP";
#use Crypt::Primes 'maurer';
use Math::Pari qw/isprime/;
use Crypt::Random 'makerandom';
-use Data::Dump 'dump';
+use Data::Dump::Filtered 'dump_filtered';
+my $bifilter = sub { my($ctx, $n) = @_;
+ return {dump=>"$n"} if ref($n) eq "Math::BigInt";
+ undef; };
$|++;
my $num = 71;
@@ -56,7 +59,7 @@ foreach my $certn (@certs) {
print proof_mark($certn->[1]);
next if $v;
print "\n\n$certn->[0] didn't verify!\n\n";
- print dump($certn->[1]);
+ print dump_filtered($certn->[1], $bifilter);
die;
}
print "\n";
@@ -64,8 +67,8 @@ 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'; }
+ if (!defined $type) { die "\nNo cert:\n\n", dump_filtered($cert, $bifilter);
}
+ if ($type =~ /n-1/i) { return ($cert->[2]->[0] eq 'B') ? '7' : '5';
}
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
