From 487de07a0a3bad1b47206a54968c2b33bde81618 Mon Sep 17 00:00:00 2001 From: Peter John Acklam Date: Sun, 2 Jan 2011 13:13:31 -0800 Subject: [perl #81400] Fix bmodinv() part of RT 63237 The following standard definition is used: z is the modular inverse of x (mod y) if and only if x*z (mod y) = 1 (mod y). - dist/Math-BigInt/lib/Math/BigInt.pm: Fix the code in bmodinv() so it can handle negative arguments. The code can be optimized further for speed, but correctnes first. - dist/Math-BigInt/t/bigintpm.inc: Fix the test case for modinv(-3, -5). The output should be -3, since -3 * -3 (mod -5) = -9 (mod -5) = -4, and 1 (mod -5) = -4. - dist/Math-BigRat/t/bigratpm.inc: Fix same test case as above. Math::BigRat::bmodinv() only handles integers, and is essentially just a front-end to Math::BigInt::bmodinv(). --- dist/Math-BigInt/lib/Math/BigInt.pm | 97 +++++++++++++++++++++++++++---------- dist/Math-BigInt/t/bigintpm.inc | 2 +- dist/Math-BigRat/t/bigratpm.inc | 2 +- 3 files changed, 74 insertions(+), 27 deletions(-) (limited to 'dist') diff --git a/dist/Math-BigInt/lib/Math/BigInt.pm b/dist/Math-BigInt/lib/Math/BigInt.pm index 3e9196a985..52acc7a682 100644 --- a/dist/Math-BigInt/lib/Math/BigInt.pm +++ b/dist/Math-BigInt/lib/Math/BigInt.pm @@ -1785,10 +1785,12 @@ sub bmod sub bmodinv { - # Modular inverse. given a number which is (hopefully) relatively - # prime to the modulus, calculate its inverse using Euclid's - # algorithm. If the number is not relatively prime to the modulus - # (i.e. their gcd is not one) then NaN is returned. + # Return modular multiplicative inverse: z is the modular inverse of x (mod + # y) if and only if x*z (mod y) = 1 (mod y). If the modulus y is larger than + # one, x and z are relative primes (i.e., their greatest common divisor is + # one). + # + # If no modular multiplicative inverse exists, NaN is returned. # set up parameters my ($self,$x,$y,@r) = (undef,@_); @@ -1800,22 +1802,58 @@ sub bmodinv return $x if $x->modify('bmodinv'); - return $x->bnan() - if ($y->{sign} ne '+' # -, NaN, +inf, -inf - || $x->is_zero() # or num == 0 - || $x->{sign} !~ /^[+-]$/ # or num NaN, inf, -inf - ); - - # put least residue into $x if $x was negative, and thus make it positive - $x->bmod($y) if $x->{sign} eq '-'; - - my $sign; - ($x->{value},$sign) = $CALC->_modinv($x->{value},$y->{value}); - return $x->bnan() if !defined $x->{value}; # in case no GCD found - return $x if !defined $sign; # already real result - $x->{sign} = $sign; # flip/flop see below - $x->bmod($y); # calc real result - $x; + # Return NaN if one or both arguments is +inf, -inf, or nan. + + return $x->bnan() if ($y->{sign} !~ /^[+-]$/ || + $x->{sign} !~ /^[+-]$/); + + # Return NaN if $y is zero; 1 % 0 makes no sense. + + return $x->bnan() if $y->is_zero(); + + # Return 0 in the trivial case. $x % 1 or $x % -1 is zero for all finite + # integers $x. + + return $x->bzero() if ($y->is_one() || + $y->is_one('-')); + + # Return NaN if $x = 0, or $x modulo $y is zero. The only valid case when + # $x = 0 is when $y = 1 or $y = -1, but that was covered above. + # + # Note that computing $x modulo $y here affects the value we'll feed to + # $CALC->_modinv() below when $x and $y have opposite signs. E.g., if $x = + # 5 and $y = 7, those two values are fed to _modinv(), but if $x = -5 and + # $y = 7, the values fed to _modinv() are $x = 2 (= -5 % 7) and $y = 7. + # The value if $x is affected only when $x and $y have opposite signs. + + $x->bmod($y); + return $x->bnan() if $x->is_zero(); + + # Compute the modular multiplicative inverse of the absolute values. We'll + # correct for the signs of $x and $y later. Return NaN if no GCD is found. + + ($x->{value}, $x->{sign}) = $CALC->_modinv($x->{value}, $y->{value}); + return $x->bnan() if !defined $x->{value}; + + # When one or both arguments are negative, we have the following + # relations. If x and y are positive: + # + # modinv(-x, -y) = -modinv(x, y) + # modinv(-x, y) = y - modinv(x, y) = -modinv(x, y) (mod y) + # modinv( x, -y) = modinv(x, y) - y = modinv(x, y) (mod -y) + + # We must swap the sign of the result if the original $x is negative. + # However, we must compensate for ignoring the signs when computing the + # inverse modulo. The net effect is that we must swap the sign of the + # result if $y is negative. + + $x -> bneg() if $y->{sign} eq '-'; + + # Compute $x modulo $y again after correcting the sign. + + $x -> bmod($y) if $x->{sign} ne $y->{sign}; + + return $x; } sub bmodpow @@ -3175,7 +3213,7 @@ Math::BigInt - Arbitrary size integer/float math package $x->bmod($y); # modulus (x % y) $x->bmodpow($exp,$mod); # modular exponentation (($num**$exp) % $mod)) - $x->bmodinv($mod); # the inverse of $x in the given modulus $mod + $x->bmodinv($mod); # the multiplicative inverse of $x modulo $mod $x->bpow($y); # power of arguments (x ** y) $x->blsft($y); # left shift in base 2 @@ -3694,11 +3732,20 @@ This method was added in v1.87 of Math::BigInt (June 2007). =head2 bmodinv() - num->bmodinv($mod); # modular inverse + $x->bmodinv($mod); # modular multiplicative inverse + +Returns the multiplicative inverse of C<$x> modulo C<$mod>. If + + $y = $x -> copy() -> bmodinv($mod) + +then C<$y> is the number closest to zero, and with the same sign as C<$mod>, +satisfying + + ($x * $y) % $mod = 1 % $mod -Returns the inverse of C<$num> in the given modulus C<$mod>. 'C' is -returned unless C<$num> is relatively prime to C<$mod>, i.e. unless -C. +If C<$x> and C<$y> are non-zero, they must be relative primes, i.e., +C. 'C' is returned when no modular multiplicative +inverse exists. =head2 bmodpow() diff --git a/dist/Math-BigInt/t/bigintpm.inc b/dist/Math-BigInt/t/bigintpm.inc index e8460041a6..822b44610e 100644 --- a/dist/Math-BigInt/t/bigintpm.inc +++ b/dist/Math-BigInt/t/bigintpm.inc @@ -1704,13 +1704,13 @@ abc:5:NaN # bmodinv Expected Results from normal use 1:5:1 3:5:2 +3:-5:-3 -2:5:2 8:5033:4404 1234567891:13:6 -1234567891:13:7 324958749843759385732954874325984357439658735983745:2348249874968739:1741662881064902 ## bmodinv Error cases / useless use of function -3:-5:NaN inf:5:NaN 5:inf:NaN -inf:5:NaN diff --git a/dist/Math-BigRat/t/bigratpm.inc b/dist/Math-BigRat/t/bigratpm.inc index 13ec6240d8..de5e9e1687 100644 --- a/dist/Math-BigRat/t/bigratpm.inc +++ b/dist/Math-BigRat/t/bigratpm.inc @@ -185,13 +185,13 @@ abc:5:NaN # bmodinv Expected Results from normal use 1:5:1 3:5:2 +3:-5:-3 -2:5:2 8:5033:4404 1234567891:13:6 -1234567891:13:7 324958749843759385732954874325984357439658735983745:2348249874968739:1741662881064902 ## bmodinv Error cases / useless use of function -3:-5:NaN inf:5:NaN 5:inf:NaN -inf:5:NaN -- cgit v1.2.1