[libmath-prime-util-perl] 121/181: Speed up some testing. Move znorder to PP.pm

Partha P. Mukherjee ppm-guest at moszumanska.debian.org
Thu May 21 18:51:13 UTC 2015


This is an automated email from the git hooks/post-receive script.

ppm-guest pushed a commit to annotated tag v0.36
in repository libmath-prime-util-perl.

commit 0c0f9d4aec20be6543aace1bbf5f62fbae7a63b9
Author: Dana Jacobsen <dana at acm.org>
Date:   Sun Jan 5 18:49:46 2014 -0800

    Speed up some testing.  Move znorder to PP.pm
---
 XS.xs                     |  2 +-
 lib/Math/Prime/Util.pm    | 47 ++++-------------------------------------
 lib/Math/Prime/Util/PP.pm | 54 +++++++++++++++++++++++++++++++++++------------
 primality.c               | 10 +++++----
 t/17-pseudoprime.t        | 11 +++++-----
 t/80-pp.t                 | 29 ++++++++++++++++---------
 6 files changed, 76 insertions(+), 77 deletions(-)

diff --git a/XS.xs b/XS.xs
index 7e791ea..28f2b26 100644
--- a/XS.xs
+++ b/XS.xs
@@ -643,7 +643,7 @@ znorder(IN SV* sva, IN SV* svn)
       XSRETURN_UV(ret);
     }
     switch (ix) {
-      case 0:  _vcallsub("_generic_znorder");  break;
+      case 0:  _vcallsub("PP::znorder");  break;
       /* TODO: Fixup public PP legendre_phi */
       case 1:
       default: _vcallsub("PP::_legendre_phi"); break;
diff --git a/lib/Math/Prime/Util.pm b/lib/Math/Prime/Util.pm
index a0d401d..ba6aaef 100644
--- a/lib/Math/Prime/Util.pm
+++ b/lib/Math/Prime/Util.pm
@@ -116,7 +116,7 @@ BEGIN {
     *carmichael_lambda = \&Math::Prime::Util::_generic_carmichael_lambda;
     *kronecker     = \&Math::Prime::Util::_generic_kronecker;
     *divisor_sum   = \&Math::Prime::Util::_generic_divisor_sum;
-    *znorder       = \&Math::Prime::Util::_generic_znorder;
+    *znorder       = \&Math::Prime::Util::PP::znorder;
     *znprimroot    = \&Math::Prime::Util::_generic_znprimroot;
     *legendre_phi  = \&Math::Prime::Util::PP::_legendre_phi;
     *gcd           = \&Math::Prime::Util::PP::gcd;
@@ -1491,46 +1491,6 @@ sub _generic_carmichael_lambda {
   return $lcm;
 }
 
-sub _generic_znorder {
-  my($a, $n) = @_;
-  _validate_num($a) || _validate_positive_integer($a);
-  _validate_num($n) || _validate_positive_integer($n);
-  return if $n <= 0;
-  return (undef,1)[$a] if $a <= 1;
-  return 1 if $n == 1;
-
-  # Sadly, Calc/FastCalc are horrendously slow for this function.
-  return if Math::BigInt::bgcd($a, $n) > 1;
-
-  # The answer is one of the divisors of phi(n) and lambda(n).
-  my $lambda = carmichael_lambda($n);
-  $a = Math::BigInt->new("$a") unless ref($a) eq 'Math::BigInt';
-
-  # This is easy and usually fast, but can bog down with too many divisors.
-  if ($lambda <= $_XS_MAXVAL) {
-    foreach my $k (divisors($lambda)) {
-      return $k if $a->copy->bmodpow("$k", $n)->is_one;
-    }
-    return;
-  }
-
-  # Algorithm 1.7 from A. Das applied to Carmichael Lambda.
-  $lambda = Math::BigInt->new("$lambda") unless ref($lambda) eq 'Math::BigInt';
-  my $k = Math::BigInt->bone;
-  foreach my $f (factor_exp($lambda)) {
-    my($pi, $ei, $enum) = (Math::BigInt->new("$f->[0]"), $f->[1], 0);
-    my $phidiv = $lambda / ($pi**$ei);
-    my $b = $a->copy->bmodpow($phidiv, $n);
-    while ($b != 1) {
-      return if $enum++ >= $ei;
-      $b = $b->copy->bmodpow($pi, $n);
-      $k *= $pi;
-    }
-  }
-  $k = _bigint_to_int($k) if $k->bacmp(''.~0) <= 0;
-  return $k;
-}
-
 sub _generic_znprimroot {
   my($n) = @_;
   $n = -$n if defined $n && $n =~ /^-\d+/;   # TODO: fix this for string bigints
@@ -1703,10 +1663,11 @@ sub lucas_sequence {
     if ref($_[0]) ne 'Math::BigInt' && $n <= $_XS_MAXVAL
     && ref($_[3]) ne 'Math::BigInt' && $k <= $_XS_MAXVAL;
   if ($_HAVE_GMP && defined &Math::Prime::Util::GMP::lucas_sequence) {
-    return map { ($_ > ~0) ? Math::BigInt->new(''.$_) : $_ }
+    return map { ($_ > ''.~0) ? Math::BigInt->new(''.$_) : $_ }
            Math::Prime::Util::GMP::lucas_sequence($n, $P, $Q, $k);
   }
-  return Math::Prime::Util::PP::lucas_sequence($n, $P, $Q, $k);
+  return map { ($_ <= ''.~0) ? _bigint_to_int($_) : $_ }
+         Math::Prime::Util::PP::lucas_sequence($n, $P, $Q, $k);
 }
 
 sub miller_rabin_random {
diff --git a/lib/Math/Prime/Util/PP.pm b/lib/Math/Prime/Util/PP.pm
index 211b2a4..a686551 100644
--- a/lib/Math/Prime/Util/PP.pm
+++ b/lib/Math/Prime/Util/PP.pm
@@ -940,18 +940,6 @@ sub _is_perfect_power {
   0;
 }
 
-sub _order {
-  my($r, $n, $lim) = @_;
-  $lim = $r unless defined $lim;
-
-  return 1 if ($n % $r) == 1;
-  for (my $j = 2; $j <= $lim; $j++) {
-    return $j if _powmod($n, $j, $r) == 1;
-  }
-  return $lim+1;
-}
-
-
 
 sub is_pseudoprime {
   my($n, $base) = @_;
@@ -1104,6 +1092,46 @@ sub _is_perfect_square {
   0;
 }
 
+sub znorder {
+  my($a, $n) = @_;
+  _validate_num($a) || _validate_positive_integer($a);
+  _validate_num($n) || _validate_positive_integer($n);
+  return if $n <= 0;
+  return (undef,1)[$a] if $a <= 1;
+  return 1 if $n == 1;
+
+  # Sadly, Calc/FastCalc are horrendously slow for this function.
+  return if Math::BigInt::bgcd($a, $n) > 1;
+
+  # The answer is one of the divisors of phi(n) and lambda(n).
+  my $lambda = Math::Prime::Util::carmichael_lambda($n);
+  $a = Math::BigInt->new("$a") unless ref($a) eq 'Math::BigInt';
+
+  # This is easy and usually fast, but can bog down with too many divisors.
+  if ($lambda <= 2**64) {
+    foreach my $k (Math::Prime::Util::divisors($lambda)) {
+      return $k if $a->copy->bmodpow("$k", $n)->is_one;
+    }
+    return;
+  }
+
+  # Algorithm 1.7 from A. Das applied to Carmichael Lambda.
+  $lambda = Math::BigInt->new("$lambda") unless ref($lambda) eq 'Math::BigInt';
+  my $k = Math::BigInt->bone;
+  foreach my $f (Math::Prime::Util::factor_exp($lambda)) {
+    my($pi, $ei, $enum) = (Math::BigInt->new("$f->[0]"), $f->[1], 0);
+    my $phidiv = $lambda / ($pi**$ei);
+    my $b = $a->copy->bmodpow($phidiv, $n);
+    while ($b != 1) {
+      return if $enum++ >= $ei;
+      $b = $b->copy->bmodpow($pi, $n);
+      $k *= $pi;
+    }
+  }
+  $k = _bigint_to_int($k) if $k->bacmp(''.~0) <= 0;
+  return $k;
+}
+
 # Find first D in sequence (5,-7,9,-11,13,-15,...) where (D|N) == -1
 sub _lucas_selfridge_params {
   my($n) = @_;
@@ -1508,7 +1536,7 @@ sub is_aks_prime {
   while ($r < $n) {
     return 0 if !($n % $r);
     #return 1 if $r >= $sqrtn;
-    last if _order($r, $n, $limit) > $limit;
+    last if znorder($r, $n) > $limit;
     $r = next_prime($r);
   }
 
diff --git a/primality.c b/primality.c
index 1a7b532..75dca8b 100644
--- a/primality.c
+++ b/primality.c
@@ -290,11 +290,12 @@ int _XS_BPSW(UV const n)
       V = montP;
       { UV v = d; b = 1; while (v >>= 1) b++; }
       while (b-- > 1) {
+        UV T = submod(  mont_prod64(V, W, n, npi),  montP, n);
         if ( (d >> (b-1)) & UVCONST(1) ) {
-          V = submod(  mont_prod64(V, W, n, npi),  montP, n);
+          V = T;
           W = submod(  mont_prod64(W, W, n, npi),  mont2, n);
         } else {
-          W = submod(  mont_prod64(V, W, n, npi),  montP, n);
+          W = T;
           V = submod(  mont_prod64(V, V, n, npi),  mont2, n);
         }
       }
@@ -628,11 +629,12 @@ int _XS_is_almost_extra_strong_lucas_pseudoprime(UV n, UV increment)
   W = mulsubmod(P, P, 2, n);
   V = P;
   while (b--) {
+    UV T = mulsubmod(V, W, P, n);
     if ( (d >> b) & UVCONST(1) ) {
-      V = mulsubmod(V, W, P, n);
+      V = T;
       W = mulsubmod(W, W, 2, n);
     } else {
-      W = mulsubmod(V, W, P, n);
+      W = T;
       V = mulsubmod(V, V, 2, n);
     }
   }
diff --git a/t/17-pseudoprime.t b/t/17-pseudoprime.t
index 311c0d1..7f6ed12 100644
--- a/t/17-pseudoprime.t
+++ b/t/17-pseudoprime.t
@@ -14,6 +14,7 @@ use Math::Prime::Util qw/is_prime
                          lucas_sequence/;
 
 my $use64 = Math::Prime::Util::prime_get_config->{'maxbits'} > 32;
+my $usexs = Math::Prime::Util::prime_get_config->{'xs'};
 my $extra = defined $ENV{EXTENDED_TESTING} && $ENV{EXTENDED_TESTING};
 
 # small primes
@@ -207,18 +208,16 @@ if ($extra) {
   is($mr2fail, 0, "is_strong_pseudoprime bases 2,3 matches is_prime");
 }
 
-# Lucas sequences, used for quite a few tests
-sub lucas_sequence_to_native {
-  map { (ref($_) eq 'Math::BigInt') ? int($_->bstr) : $_ } lucas_sequence(@_);
-}
+# Lucas sequences, used for quite a few primality tests
 while (my($params, $expect) = each (%lucas_sequences)) {
-  is_deeply( [lucas_sequence_to_native(split(' ', $params))], $expect, "Lucas sequence $params" );
+  is_deeply( [lucas_sequence(split(' ', $params))], $expect, "Lucas sequence $params" );
 }
 
 SKIP: {
   skip "Old Perl+bigint segfaults in F-U code",1+2*$use64 if $] < 5.008;
   my $fufail = 0;
-  foreach my $i (1 .. 100) {
+  my $ntests = ($usexs) ? 100 : 2;
+  foreach my $i (1 .. $ntests) {
     my $n = 2*int(rand(1000000000)) + 1;
     my $ispfu = !!is_frobenius_underwood_pseudoprime($n);
     my $prime = !!is_prime($n);
diff --git a/t/80-pp.t b/t/80-pp.t
index c8b0c59..bf822e3 100644
--- a/t/80-pp.t
+++ b/t/80-pp.t
@@ -254,13 +254,15 @@ plan tests => 2 +
               4 + scalar(keys %pseudoprimes) +
               scalar(keys %eivals) + scalar(keys %livals) + scalar(keys %rvals) + scalar(keys %rzvals) +
               ($extra ? 4 : 0) +  # Bigfloat RiemannZeta
-              1 + 1 +    # factor
-              10 + 8*3 + # factoring subs
-              10 +       # AKS
+              1 + 1 +             # factor
+              10 + 7*3 +          # factoring subs
+              1 +                 # HOLF
+              ($extra ? 3 : 0) +  # HOLF extra
+              10 +                # AKS
               ($use64 ? 3 : 2) +  # Lucas and BLS75 primality proofs
-              4 +        # M-R and Lucas on bigint
-              13 +       # Misc util.pm functions
-              scalar(keys %ipp) +
+              4 +                 # M-R and Lucas on bigint
+              13 +                # Misc util.pm functions
+              scalar(keys %ipp) + # is_prob_prime
               1;
 
 use Math::Prime::Util qw/primes prime_count_approx prime_count_lower
@@ -529,10 +531,6 @@ if ($extra) {
   # 1013 4294967197 4294967291
   my $nbig = Math::BigInt->new("18686551294184381720251");
   my @nfac;
-  @nfac = sort {$a<=>$b} Math::Prime::Util::PP::holf_factor($nbig);
-  is(scalar @nfac, 2, "holf finds a factor of 18686551294184381720251");
-  is($nfac[0] * $nfac[1], $nbig, "holf found a correct factor");
-  ok($nfac[0] != 1 && $nfac[1] != 1, "holf didn't return a degenerate factor");
   @nfac = sort {$a<=>$b} Math::Prime::Util::PP::prho_factor($nbig);
   is(scalar @nfac, 2, "prho finds a factor of 18686551294184381720251");
   is($nfac[0] * $nfac[1], $nbig, "prho found a correct factor");
@@ -563,6 +561,17 @@ if ($extra) {
   is(scalar @nfac, 2, "fermat finds a factor of 73786976930493367637");
   is($nfac[0] * $nfac[1], $nbig, "fermat found a correct factor");
   ok($nfac[0] != 1 && $nfac[1] != 1, "fermat didn't return a degenerate factor");
+  if ($extra) {
+    @nfac = sort {$a<=>$b} Math::Prime::Util::PP::holf_factor($nbig);
+    is(scalar @nfac, 2, "holf finds a factor of 18686551294184381720251");
+    is($nfac[0] * $nfac[1], $nbig, "holf found a correct factor");
+    ok($nfac[0] != 1 && $nfac[1] != 1, "holf didn't return a degenerate factor");
+  }
+  {
+    $nbig = Math::BigInt->new("99999999999979999998975857");
+    @nfac = sort {$a<=>$b} Math::Prime::Util::PP::holf_factor($nbig);
+    is_deeply(\@nfac, [9999999998987,10000000001011], "holf correctly factors 99999999999979999998975857");
+  }
   SKIP: {
     # Unfortunately we can't guarantee this isn't found in stage 1.
     skip "ecm stage 2", 3 unless $extra;

-- 
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-perl/packages/libmath-prime-util-perl.git



More information about the Pkg-perl-cvs-commits mailing list