[libmath-prime-util-perl] 54/181: Huge speedup for znprimroot with non-cyclic input

[Partha P. Mukherjee]

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 7033eec7f3f75bca50658ce6846e7ce18ee31828
Author: Dana Jacobsen <dana at acm.org>
Date:   Sat Dec 28 12:16:53 2013 -0800

    Huge speedup for znprimroot with non-cyclic input
---
 lib/Math/Prime/Util.pm | 25 +++++++++++++++++++------
 t/19-moebius.t         |  5 ++++-
 util.c                 | 10 +++++++---
 3 files changed, 30 insertions(+), 10 deletions(-)

diff --git a/lib/Math/Prime/Util.pm b/lib/Math/Prime/Util.pm
index e3d9468..64de45b 100644
--- a/lib/Math/Prime/Util.pm
+++ b/lib/Math/Prime/Util.pm
@@ -1830,8 +1830,8 @@ sub _generic_kronecker {
   my($a, $b) = @_;
   croak "Parameter must be defined" if !defined $a;
   croak "Parameter must be defined" if !defined $b;
-  croak "Parameter '$a' must be an integer" unless $a =~ /^-?\d+/;
-  croak "Parameter '$b' must be an integer" unless $b =~ /^-?\d+/;
+  croak "Parameter '$a' must be an integer" unless $a =~ /^[-+]?\d+/;
+  croak "Parameter '$b' must be an integer" unless $b =~ /^[-+]?\d+/;
 
   return Math::BigInt->new(''.Math::Prime::Util::GMP::kronecker($a,$b))
     if $_HAVE_GMP && defined &Math::Prime::Util::GMP::kronecker;
@@ -3860,9 +3860,6 @@ the primitive root exists, while L<OEIS A046145|http://oeis.org/A046145>
 is a list of the smallest primitive roots, which is what this function
 produces.
 
-This will always produce the smallest primitive root.  Pari does not
-necessarily produce the smallest result for composites.
-
 
 =head2 random_prime
 
@@ -4556,6 +4553,17 @@ Compute L<OEIS A054903|http://oeis.org/A054903> just like CRG4's Pari example:
     say if divisor_sum($_)+6 == divisor_sum($_+6)
   } 9,1e7;
 
+Construct the table shown in L<OEIS A046147|http://oeis.org/A046147>:
+
+  use Math::Prime::Util qw/znorder euler_phi/;
+  foreach my $n (1..100) {
+    my $phi = euler_phi($n);
+    my @r = grep {    Math::BigInt::bgcd($_,$n) == 1
+                   && znorder($_,$n) == $phi
+                 } 1..$n-1;
+    say "$n ", join(" ", @r);
+  }
+
 =head1 PRIMALITY TESTING NOTES
 
 Above C<2^64>, L</is_prob_prime> performs an extra-strong
@@ -4853,7 +4861,12 @@ Math::Pari.
 
 =item C<kronecker>, C<znorder>, C<znprimroot>
 
-Similar to MPU's L</kronecker>, L</znorder>, and L</znprimroot>.
+Similar to MPU's L</kronecker>, L</znorder>, and L</znprimroot>.  Pari's
+C<znprimroot> only returns the smallest root for prime powers.  The
+behavior is undefined when the group is not cyclic (sometimes it throws
+an exception, sometimes it returns an incorrect answer).
+MPU's L</znprimroot> will always return the smallest root if it exists,
+and C<undef> otherwise.
 
 =item C<sigma>
 
diff --git a/t/19-moebius.t b/t/19-moebius.t
index 1423ccb..328846f 100644
--- a/t/19-moebius.t
+++ b/t/19-moebius.t
@@ -214,6 +214,7 @@ my %primroots = (
 1407827621 =>   2,
 1520874431 =>  17,
 1685283601 => 164,
+ 100000001 => undef,
 );
 if ($use64) {
   $primroots{2232881419280027} = 6;
@@ -253,7 +254,7 @@ plan tests => 0 + 1
                 + 7 + scalar(keys %totients)
                 + 1 # Small Carmichael Lambda
                 + scalar(@mult_orders)
-                + scalar(keys %primroots)
+                + scalar(keys %primroots) + 2
                 + scalar(keys %jordan_totients)
                 + 2  # Dedekind psi calculated two ways
                 + 2  # Calculate J5 two different ways
@@ -404,6 +405,8 @@ foreach my $moarg (@mult_orders) {
 while (my($n, $root) = each (%primroots)) {
   is( znprimroot($n), $root, "znprimroot($n) == " . ((defined $root) ? $root : "<undef>") );
 }
+is( znprimroot("-100000898"), 31, "znprimroot(\"-100000898\") == 31" );
+is( znprimroot("+100000898"), 31, "znprimroot(\"+100000898\") == 31" );
 ###### liouville
 foreach my $i (@liouville_pos) {
   is( liouville($i),  1, "liouville($i) = 1" );
diff --git a/util.c b/util.c
index d2c18e7..336a333 100644
--- a/util.c
+++ b/util.c
@@ -1015,17 +1015,21 @@ UV znprimroot(UV n) {
   UV a, phi;
   int i, nfactors;
   if (n <= 4) return (n == 0) ? 0 : n-1;
+  if (n % 4 == 0)  return 0;
   phi = totient(n);
+  /* Check if a primitive root exists. */
+  if (!_XS_is_prob_prime(n) && phi != carmichael_lambda(n))  return 0;
   nfactors = factor_exp(phi, fac, exp);
   for (i = 0; i < nfactors; i++)
     exp[i] = phi / fac[i];  /* exp[i] = phi(n) / i-th-factor-of-phi(n) */
   for (a = 2; a < n; a++) {
     if (kronecker_uu(a, n) == 0)  continue;
-    for (i = 0; i < nfactors; i++) {
-      if (powmod(a, exp[i], n) == 1) break;
-    }
+    for (i = 0; i < nfactors; i++)
+      if (powmod(a, exp[i], n) == 1)
+        break;
     if (i == nfactors) return a;
   }
+  return 0;
 }
 
 

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