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

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

commit 174b33419add2a4e9162f4db5969741d2777950f
Author: Dana Jacobsen <d...@acm.org>
Date:   Mon May 5 12:36:43 2014 -0700

    invmod edge cases, add tests
---
 TODO                        |  2 --
 XS.xs                       | 67 ++++++++++++++++++++++++++-------------------
 lib/Math/Prime/Util.pm      | 13 ++++++---
 lib/Math/Prime/Util/PP.pm   |  5 +++-
 lib/Math/Prime/Util/PPFE.pm |  6 ++++
 t/19-moebius.t              | 28 +++++++++++++++++++
 t/92-release-pod-coverage.t |  2 +-
 7 files changed, 87 insertions(+), 36 deletions(-)

diff --git a/TODO b/TODO
index af89a74..f8132c7 100644
--- a/TODO
+++ b/TODO
@@ -72,5 +72,3 @@
 - We don't use legendre_phi for other functions any more, but it'd be nice
   to speed it up using some ideas from the Ohana 2011 SAGE branch.  For example
   (10**13,10**5) takes 2.5x longer, albeit with 6x less memory.
-
-- tests for invmod, inc negatives, undef, non-numeric, high-bit, and no inv.
diff --git a/XS.xs b/XS.xs
index 5dba66b..0e21827 100644
--- a/XS.xs
+++ b/XS.xs
@@ -737,9 +737,8 @@ divisor_sum(IN SV* svn, ...)
 void
 znorder(IN SV* sva, IN SV* svn)
   ALIAS:
-    invmod = 1
-    jordan_totient = 2
-    legendre_phi = 3
+    jordan_totient = 1
+    legendre_phi = 2
   PREINIT:
     int astatus, nstatus;
   PPCODE:
@@ -752,25 +751,22 @@ znorder(IN SV* sva, IN SV* svn)
       switch (ix) {
         case 0:  ret = znorder(a, n);
                  break;
-        case 1:  ret = modinverse(a, n);
-                 break;
-        case 2:  ret = jordan_totient(a, n);
+        case 1:  ret = jordan_totient(a, n);
                  if (ret == 0 && n > 1)
                    goto overflow;
                  break;
-        case 3:
+        case 2:
         default: ret = legendre_phi(a, n);
                  break;
       }
-      if (ret == 0 && (ix == 0 || ix == 1))  XSRETURN_UNDEF;  /* not defined */
+      if (ret == 0 && ix == 0)  XSRETURN_UNDEF;  /* not defined */
       XSRETURN_UV(ret);
     }
     overflow:
     switch (ix) {
       case 0:  _vcallsub_with_gmp("znorder");  break;
-      case 1:  _vcallsub_with_gmp("invmod");  break;
-      case 2:  _vcallsub_with_pp("jordan_totient");  break;
-      case 3:
+      case 1:  _vcallsub_with_pp("jordan_totient");  break;
+      case 2:
       default: _vcallsub_with_pp("legendre_phi"); break;
     }
     return; /* skip implicit PUTBACK */
@@ -796,35 +792,50 @@ void
 kronecker(IN SV* sva, IN SV* svb)
   ALIAS:
     valuation = 1
+    invmod = 2
   PREINIT:
     int astatus, bstatus, abpositive, abnegative;
   PPCODE:
     astatus = _validate_int(aTHX_ sva, 2);
     bstatus = _validate_int(aTHX_ svb, 2);
-    if (ix == 0) {
-      /* Are both a and b positive? */
-      abpositive = astatus == 1 && bstatus == 1;
-      /* Will both fit in IVs?  We should use a bitmask return. */
-      abnegative = !abpositive
-                   && (astatus != 0 && SvIOK(sva) && !SvIsUV(sva))
-                   && (bstatus != 0 && SvIOK(svb) && !SvIsUV(svb));
-      if (abpositive || abnegative) {
-        UV a = my_svuv(sva);
-        UV b = my_svuv(svb);
-        int k = (abpositive) ? kronecker_uu(a,b) : kronecker_ss(a,b);
-        RETURN_NPARITY(k);
-      }
-    } else {
-      if (astatus != 0 && bstatus != 0) {
+    if (astatus != 0 && bstatus != 0) {
+      if (ix == 0) {
+        /* Are both a and b positive? */
+        abpositive = astatus == 1 && bstatus == 1;
+        /* Will both fit in IVs?  We should use a bitmask return. */
+        abnegative = !abpositive
+                     && (SvIOK(sva) && !SvIsUV(sva))
+                     && (SvIOK(svb) && !SvIsUV(svb));
+        if (abpositive || abnegative) {
+          UV a = my_svuv(sva);
+          UV b = my_svuv(svb);
+          int k = (abpositive) ? kronecker_uu(a,b) : kronecker_ss(a,b);
+          RETURN_NPARITY(k);
+        }
+      } else if (ix == 1) {
         UV n = (astatus == -1) ? (UV)(-(my_sviv(sva))) : my_svuv(sva);
         UV k = (bstatus == -1) ? (UV)(-(my_sviv(svb))) : my_svuv(svb);
         XSRETURN_UV( valuation(n, k) );
+      } else {
+        UV a, n, ret = 0;
+        n = (bstatus != -1) ? my_svuv(svb) : (UV)(-(my_sviv(svb)));
+        if (n > 0) {
+          a = (astatus != -1) ? my_svuv(sva)
+                              : n * ((UV)(-my_sviv(sva))/n + 1) + my_sviv(sva);
+          if (a > 0) {
+            if (n == 1) XSRETURN_UV(0);
+            ret = modinverse(a, n);
+          }
+        }
+        if (ret == 0) XSRETURN_UNDEF;
+        XSRETURN_UV(ret);
       }
     }
     switch (ix) {
       case 0:  _vcallsub_with_gmp("kronecker");  break;
-      case 1:
-      default: _vcallsub_with_gmp("valuation"); break;
+      case 1:  _vcallsub_with_gmp("valuation"); break;
+      case 2:
+      default: _vcallsub_with_gmp("invmod"); break;
     }
     return; /* skip implicit PUTBACK */
 
diff --git a/lib/Math/Prime/Util.pm b/lib/Math/Prime/Util.pm
index 9306bb1..33f08fb 100644
--- a/lib/Math/Prime/Util.pm
+++ b/lib/Math/Prime/Util.pm
@@ -813,7 +813,7 @@ __END__
 
 =encoding utf8
 
-=for stopwords forprimes forcomposites fordivisors Möbius Deléglise totient 
moebius mertens liouville znorder irand primesieve uniqued k-tuples von SoE 
pari yafu fonction qui compte le nombre nombres voor PhD superset sqrt(N) 
gcd(A^M k-th (10001st primegen libtommath kronecker znprimroot znlog gcd lcm
+=for stopwords forprimes forcomposites fordivisors Möbius Deléglise totient 
moebius mertens liouville znorder irand primesieve uniqued k-tuples von SoE 
pari yafu fonction qui compte le nombre nombres voor PhD superset sqrt(N) 
gcd(A^M k-th (10001st primegen libtommath kronecker znprimroot znlog gcd lcm 
invmod
 
 
 =head1 NAME
@@ -1878,13 +1878,18 @@ follow the semantics of Mathematica, Pari, and Perl 6, 
re:
   lcm(0, n) = 0              Any zero in list results in zero return
   lcm(n,-m) = lcm(n, m)      We use the absolute values
 
-=head invmod
+=head2 invmod
 
   say "The inverse of 42 mod 2017 = ", invmod(42,2017);
 
 Given two integers C<a> and C<n>, return the inverse of C<a> modulo C<n>.
-If not defined, undef is returned.  If defined, then multiplying the return
-value by C<a>, modulo C<n>, will equal 1.
+If not defined, undef is returned.  If defined, then the return value
+multiplied by C<a> equals C<1> modulo C<n>.
+
+This results correspond to the Pari result of C<lift(Mod(1/a,n))>.  The
+semantics with respect to negative arguments match Pari.  Notably, a
+negative C<n> is negated, which is different from Math::BigInt, but in both
+cases the return value is still congruent to C<1> modulo C<n> as expected.
 
 =head2 valuation
 
diff --git a/lib/Math/Prime/Util/PP.pm b/lib/Math/Prime/Util/PP.pm
index 138ca09..023d023 100644
--- a/lib/Math/Prime/Util/PP.pm
+++ b/lib/Math/Prime/Util/PP.pm
@@ -1486,8 +1486,11 @@ sub lcm {
 }
 sub invmod {
   my($a,$n) = @_;
+  return if $n == 0 || $a == 0;
+  return 0 if $n == 1;
+  $n = -$n if $n < 0;  # Pari semantics
   if ($n > (~0>>1)) {
-    my $invmod = Math::BigInt->new("$a")->bmodinv($n);
+    my $invmod = Math::BigInt->new("$a")->bmodinv("$n");
     return if !defined $invmod || $invmod->is_nan;
     $invmod = _bigint_to_int($invmod) if $invmod->bacmp(''.~0) <= 0;
     return $invmod;
diff --git a/lib/Math/Prime/Util/PPFE.pm b/lib/Math/Prime/Util/PPFE.pm
index e631953..6927f80 100644
--- a/lib/Math/Prime/Util/PPFE.pm
+++ b/lib/Math/Prime/Util/PPFE.pm
@@ -322,6 +322,12 @@ sub lcm {
   return Math::Prime::Util::PP::lcm(@_);
 }
 sub invmod {
+  my ($a, $n) = @_;
+  my ($va, $vn) = ($a, $n);
+  $va = -$va if defined $va && $va < 0;
+  $vn = -$vn if defined $vn && $vn < 0;
+  _validate_positive_integer($va);
+  _validate_positive_integer($vn);
   return Math::Prime::Util::PP::invmod(@_);
 }
 
diff --git a/t/19-moebius.t b/t/19-moebius.t
index f9011bb..85a9432 100644
--- a/t/19-moebius.t
+++ b/t/19-moebius.t
@@ -7,6 +7,7 @@ use Math::Prime::Util
    qw/moebius mertens euler_phi jordan_totient divisor_sum exp_mangoldt
       chebyshev_theta chebyshev_psi carmichael_lambda znorder liouville
       znprimroot znlog kronecker legendre_phi gcd lcm is_power valuation
+      invmod
      /;
 
 my $extra = defined $ENV{EXTENDED_TESTING} && $ENV{EXTENDED_TESTING};
@@ -366,6 +367,23 @@ if ($use64) {
   push @{$powers{9}}, 118587876497;
 }
 
+my @invmods = (
+ [ 0, 0, undef],
+ [ 1, 0, undef],
+ [ 0, 1, undef],
+ [ 1, 1, 0],
+ [ 45, 59, 21],
+ [  42,  2017, 1969],
+ [  42, -2017, 1969],
+ [ -42,  2017, 48],
+ [ -42, -2017, 48],
+ [ 14, 28474, undef],
+);
+if ($use64) {
+ push @invmods, [ 13, 9223372036854775808, 5675921253449092805 ];
+ push @invmods, [ 14, 18446744073709551615, 17129119497016012214 ];
+}
+
 # These are slow with XS, and *really* slow with PP.
 if (!$usexs) {
   %big_mertens = map { $_ => $big_mertens{$_} }
@@ -399,6 +417,7 @@ plan tests => 0 + 1
                 + scalar(@znlogs)
                 + scalar(@legendre_sums)
                 + scalar(@valuations)
+                + 3 + scalar(@invmods)
                 + scalar(keys %powers)
                 + scalar(keys %primroots) + 2
                 + scalar(keys %jordan_totients)
@@ -609,6 +628,15 @@ foreach my $r (@valuations) {
   my($n, $k, $exp) = @$r;
   is( valuation($n, $k), $exp, "valuation($n,$k) = $exp" );
 }
+###### invmod
+ok(!eval { invmod(undef,11); }, "invmod(undef,11)");
+ok(!eval { invmod(11,undef); }, "invmod(11,undef)");
+ok(!eval { invmod('nan',11); }, "invmod('nan',11)");
+
+foreach my $r (@invmods) {
+  my($a, $n, $exp) = @$r;
+  is( invmod($a,$n), $exp, "invmod($a,$n) = ".((defined $exp)?$exp:"<undef>") 
);
+}
 
 sub cmp_closeto {
   my $got = shift;
diff --git a/t/92-release-pod-coverage.t b/t/92-release-pod-coverage.t
index eecf464..74eb36b 100644
--- a/t/92-release-pod-coverage.t
+++ b/t/92-release-pod-coverage.t
@@ -65,7 +65,7 @@ sub mpu_public_regex {
       random_maurer_prime random_maurer_prime_with_cert
       random_shawe_taylor_prime random_shawe_taylor_prime_with_cert
       primorial pn_primorial consecutive_integer_lcm
-      gcd lcm factor factor_exp all_factors divisors valuation
+      gcd lcm factor factor_exp all_factors divisors valuation invmod
       moebius mertens euler_phi jordan_totient exp_mangoldt liouville
       partitions
       chebyshev_theta chebyshev_psi

-- 
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

Reply via email to