Add gcd & abs (#75)

Author: konsumlamm

src/bigints.nim

@@ -71,6 +71,14 @@ func isZero(a: BigInt): bool {.inline.} =
       return false
   return true
 
+func abs*(a: BigInt): BigInt =
+  # Returns the absolute value of `a`.
+  runnableExamples:
+    assert abs(42.initBigInt) == 42.initBigInt
+    assert abs(-12.initBigInt) == 12.initBigInt
+  result = a
+  result.isNegative = false
+
 func unsignedCmp(a: BigInt, b: uint32): int64 =
   # ignores the sign of `a`
   # `a` and `b` are assumed to not be zero
@@ -751,6 +759,43 @@ func `divmod`*(a, b: BigInt): tuple[q, r: BigInt] =
     assert divmod(a, b) == (3.initBigInt, 2.initBigInt)
   division(result.q, result.r, a, b)
 
+func countTrailingZeroBits(a: BigInt): int =
+  var count = 0
+  for x in a.limbs:
+    if x == 0:
+      count += 32
+    else:
+      return count + countTrailingZeroBits(x)
+  return count
+
+func gcd*(a, b: BigInt): BigInt =
+  ## Returns the greatest common divisor (GCD) of `a` and `b`.
+  runnableExamples:
+    assert gcd(54.initBigInt, 24.initBigInt) == 6.initBigInt
+
+  # binary GCD algorithm
+  var
+    u = abs(a)
+    v = abs(b)
+  if u.isZero:
+    return v
+  elif v.isZero:
+    return u
+  let
+    i = countTrailingZeroBits(u)
+    j = countTrailingZeroBits(v)
+    k = min(i, j)
+  u = u shr i
+  v = v shr j
+  while true:
+    # u and v are odd
+    if u > v:
+      swap(u, v)
+    v -= u
+    if v.isZero:
+      return u shl k
+    v = v shr countTrailingZeroBits(v)
+
 
 func toSignedInt*[T: SomeSignedInt](x: BigInt): Option[T] =
   ## Converts a `BigInt` number to signed integer, if possible.

tests/tbigints.nim

@@ -461,6 +461,35 @@ proc main() =
     doAssertRaises(DivByZeroDefect): discard one mod zero
     doAssertRaises(DivByZeroDefect): discard divmod(one, zero)
 
+  block: # gcd
+    let a = "866506".initBigInt
+    let b = "140640".initBigInt
+    let two = 2.initBigInt
+    doAssert gcd(a, b) == two
+    # gcd(a, b) = gcd(b, a)
+    doAssert gcd(b, a) == two
+    # gcd(a, -b) = gcd(a, |b|)
+    doAssert gcd(-a, b) == two
+    doAssert gcd(a, -b) == two
+    doAssert gcd(-a, -b) == two
+
+    block: # properties of gcd
+      let a = "668403".initBigInt
+      let b = "753160".initBigInt
+      let c = "249115".initBigInt
+      doAssert gcd(a, b) == gcd(b, a)
+      doAssert gcd(a, zero) == a
+      doAssert gcd(a, a) == a
+      doAssert gcd(c * a, c * b) == c * gcd(a,b)
+      doAssert gcd(a, gcd(b, c)) == gcd(gcd(a, b), c)
+      doAssert gcd(a, b) == gcd(b, a mod b)
+
+  block: # abs
+    doAssert abs(zero) == zero
+    doAssert abs(-zero) == zero
+    doAssert abs(one) == one
+    doAssert abs(-one) == one
+
   block: # toSignedInt
     let
       a = initBigInt(7)
@@ -531,5 +560,6 @@ proc main() =
     doAssert pred(a, 3) == initBigInt(4)
     doAssert succ(a, 3) == initBigInt(10)
 
+
 static: main()
 main()