Changeset 3874

Timestamp:

08/11/08 04:21:33 (4 months ago)

Author:

Don Clugston

Message:

bignum: Major speedup for asymmetric multiplication: now ensures Karatsuba is used ONLY when it is actually beneficial; otherwise uses Schoolbook for such iterations. This algorithm is significantly better than anything else I've seen reported elsewhere.

Files:

• trunk/tango/math/internal/BiguintCore.d (modified) (8 diffs)

trunk/tango/math/internal/BiguintCore.d

@@ -184 +184 @@
-        *negative = true;
-        large = y; small = x;
-//        return biguintSubDifferentLength(y, x);
-    } else {
-        *negative = false;
-        large = x; small = y;
-  //      return biguintSubDifferentLength(x, y);
-    }
-    // result.length will be equal to larger length, or could decrease by 1.
@@ -240 +238 @@
-    } else return result[0..$-1];
-}
-
+
-/** General unsigned multiply routine for bigints.
- *  Sets result = x*y.
@@ -285 +283 @@
-    
-    uint half = (x.length >> 1) + (x.length & 1);
-y.length <= half
+2*y.length*y.length <= x.length*x.length
-        // UNBALANCED MULTIPLY
-
-
-
-
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-        
-
-
-
-
-
-
-
-
-
-
-
-
+
+
+
+
+
+
+
+
+
+
+
+
+
-        while (done < x.length) {
-
-
+
+
+
+
-            result[done + y.length .. done + y.length + chunksize] 
-                = partial[y.length .. y.length + chunksize];
-            simpleAddAssign(result[done .. done + y.length+chunksize], partial[0 .. y.length]);
-y.length
+chunksize
-        }
-        delete scratchbuff;
@@ -351 +384 @@
-in {    
-    assert(result.length == left.length + right.length);
-
+//
-}
-body {
@@ -378 +411 @@
-}
-
+uint subSimple(uint [] result, uint [] left, uint [] right)
+in {
+    assert(result.length == left.length);
+    assert(left.length >= right.length);
+    assert(right.length>0);
+}
+body {
+    uint carry = multibyteAddSub!('-')(result[0..right.length],
+            left[0..right.length], right, 0);
+    if (right.length < left.length) {
+        result[right.length..left.length] = left[right.length .. $];            
+        carry = multibyteIncrementAssign!('-')(result[right.length..$], carry);
+    } //else if (result.length==left.length+1) { result[$-1] = carry; carry=0; }
+    return carry;
+}
+
+
-/*  result must be larger than right.
-*/
@@ -411 +461 @@
-}
-
-//uint lastleftover = 0;
-/* Sets result = x*y, using Karatsuba multiplication.
-Valid only for balanced multiplies, and x must be longer than
-ie 2*y.length > x.length >= y.length
+x must be longer or equal to
+Valid only for balanced multiplies, where x is not shorter than y
-* It is efficient only if sqrt(2)*y.length > x.length >= y.length
-* Params:
@@ -425 +474 @@
-    if (x.length <= KARATSUBALIMIT) {
-        return mulSimple(result, x, y);
-
+
+
+
+
-    // Karatsuba multiply uses the following result:
-    // (Nx1 + x0)*(Ny1 + y0) = (N*N) x1y1 + x0y0 + N * mid
-    // where mid = (x1+x0)*(y1+y0) - x1y1 - x0y0
-    // requiring 3 multiplies of length N, instead of 4.
+    
+    // TODO: Knuth's variant is superior:
+    // (Nx1 + x0)*(Ny1 + y0) = (N*N) x1y1 + x0y0 - N * mid
+    // where mid = (x0-x1)*(y0-y1) - x1y1 - x0y0
+    // since then the lengths cannot exceed length N.
-
-    // half length, round up.
-    uint half = (x.length >> 1) + (x.length & 1);
-    assert(y.length>half, "Asymmetric Karatsuba");
-    
-    uint [] x0 = x[0 .. half];
@@ -460 +516 @@
-    mulKaratsuba(resultLow, x0, y0, newscratchbuff);
-    mid.simpleSubAssign(resultLow);
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-    mid.simpleSubAssign(resultHigh);
-
+    
-    // result += MID
-    result[half..$].simpleAddAssign(mid);