First step to fix PR2088. Implement routine to compute the

multiplicative inverse of a given number. Modify udivrem to allow input and 
output pairs of arguments to overlap. Patch is based on the work by Chandler
Carruth.


git-svn-id: https://llvm.org/svn/llvm-project/llvm/trunk@52638 91177308-0d34-0410-b5e6-96231b3b80d8

1 files changed, 48 insertions, 7 deletions

diff --git a/lib/Support/APInt.cpp b/lib/Support/APInt.cpp
index 38b379005a..b0577b3680 100644
--- a/lib/Support/APInt.cpp
+++ b/lib/Support/APInt.cpp
@@ -1426,6 +1426,50 @@ APInt APInt::sqrt() const {
   return x_old + 1;
 }
 
+/// Computes the multiplicative inverse of this APInt for a given modulo. The
+/// iterative extended Euclidean algorithm is used to solve for this value,
+/// however we simplify it to speed up calculating only the inverse, and take
+/// advantage of div+rem calculations. We also use some tricks to avoid copying
+/// (potentially large) APInts around.
+APInt APInt::multiplicativeInverse(const APInt& modulo) const {
+  assert(ult(modulo) && "This APInt must be smaller than the modulo");
+
+  // Using the properties listed at the following web page (accessed 06/21/08):
+  //   http://www.numbertheory.org/php/euclid.html
+  // (especially the properties numbered 3, 4 and 9) it can be proved that
+  // BitWidth bits suffice for all the computations in the algorithm implemented
+  // below. More precisely, this number of bits suffice if the multiplicative
+  // inverse exists, but may not suffice for the general extended Euclidean
+  // algorithm.
+
+  APInt r[2] = { modulo, *this };
+  APInt t[2] = { APInt(BitWidth, 0), APInt(BitWidth, 1) };
+  APInt q(BitWidth, 0);
+  
+  unsigned i;
+  for (i = 0; r[i^1] != 0; i ^= 1) {
+    // An overview of the math without the confusing bit-flipping:
+    // q = r[i-2] / r[i-1]
+    // r[i] = r[i-2] % r[i-1]
+    // t[i] = t[i-2] - t[i-1] * q
+    udivrem(r[i], r[i^1], q, r[i]);
+    t[i] -= t[i^1] * q;
+  }
+
+  // If this APInt and the modulo are not coprime, there is no multiplicative
+  // inverse, so return 0. We check this by looking at the next-to-last
+  // remainder, which is the gcd(*this,modulo) as calculated by the Euclidean
+  // algorithm.
+  if (r[i] != 1)
+    return APInt(BitWidth, 0);
+
+  // The next-to-last t is the multiplicative inverse.  However, we are
+  // interested in a positive inverse. Calcuate a positive one from a negative
+  // one if necessary. A simple addition of the modulo suffices because
+  // abs(t[i]) is known to less than *this/2 (see the link above).
+  return t[i].isNegative() ? t[i] + modulo : t[i];
+}
+
 /// Implementation of Knuth's Algorithm D (Division of nonnegative integers)
 /// from "Art of Computer Programming, Volume 2", section 4.3.1, p. 272. The
 /// variables here have the same names as in the algorithm. Comments explain
@@ -1882,13 +1926,10 @@ void APInt::udivrem(const APInt &LHS, const APInt &RHS,
   
   if (lhsWords == 1 && rhsWords == 1) {
     // There is only one word to consider so use the native versions.
-    if (LHS.isSingleWord()) {
-      Quotient = APInt(LHS.getBitWidth(), LHS.VAL / RHS.VAL);
-      Remainder = APInt(LHS.getBitWidth(), LHS.VAL % RHS.VAL);
-    } else {
-      Quotient = APInt(LHS.getBitWidth(), LHS.pVal[0] / RHS.pVal[0]);
-      Remainder = APInt(LHS.getBitWidth(), LHS.pVal[0] % RHS.pVal[0]);
-    }
+    uint64_t lhsValue = LHS.isSingleWord() ? LHS.VAL : LHS.pVal[0];
+    uint64_t rhsValue = RHS.isSingleWord() ? RHS.VAL : RHS.pVal[0];
+    Quotient = APInt(LHS.getBitWidth(), lhsValue / rhsValue);
+    Remainder = APInt(LHS.getBitWidth(), lhsValue % rhsValue);
     return;
   }