Add exponentiation - closes #37 (#46)

* Add exponentiation

* Change name to pow (like ttmath) - `^` to discuss

stint/private/uint_exp.nim

@@ -0,0 +1,31 @@
+# Stint
+# Copyright 2018 Status Research & Development GmbH
+# Licensed under either of
+#
+#  * Apache License, version 2.0, ([LICENSE-APACHE](LICENSE-APACHE) or http://www.apache.org/licenses/LICENSE-2.0)
+#  * MIT license ([LICENSE-MIT](LICENSE-MIT) or http://opensource.org/licenses/MIT)
+#
+# at your option. This file may not be copied, modified, or distributed except according to those terms.
+
+import
+  ./datatypes,
+  ./uint_bitwise_ops, ./uint_mul, ./initialization
+
+func pow*(x: UintImpl, y: Natural): UintImpl =
+  ## Compute ``x`` to the power of ``y``,
+  ## ``x`` must be non-negative
+
+  # Implementation uses exponentiation by squaring
+  # See Nim math module: https://github.com/nim-lang/Nim/blob/4ed24aa3eb78ba4ff55aac3008ec3c2427776e50/lib/pure/math.nim#L429
+  # And Eli Bendersky's blog: https://eli.thegreenplace.net/2009/03/21/efficient-integer-exponentiation-algorithms
+
+  var (x, y) = (x, y)
+  result = one(type x)
+
+  while true:
+    if (y and 1) != 0:
+      result = result * x
+    y = y shr 1
+    if y == 0:
+      break
+    x = x * x

stint/uint_public.nim

@@ -61,9 +61,9 @@ make_unary(`not`, StUint)
 make_binary(`or`, StUint)
 make_binary(`and`, StUint)
 make_binary(`xor`, StUint)
-proc `shr`*(x: StUint, y: SomeInteger): StUint {.inline, noSideEffect.} =
+func `shr`*(x: StUint, y: SomeInteger): StUint {.inline.} =
   result.data = x.data shr y
-proc `shl`*(x: StUint, y: SomeInteger): StUint {.inline, noSideEffect.} =
+func `shl`*(x: StUint, y: SomeInteger): StUint {.inline.} =
   result.data = x.data shl y
 
 import ./private/uint_highlow
@@ -86,3 +86,11 @@ func one*[bits: static[int]](T: typedesc[Stuint[bits] or Stint[bits]]): T {.inli
 
 func zero*[bits: static[int]](T: typedesc[Stuint[bits] or Stint[bits]]): T {.inline.} =
   discard
+
+import ./private/uint_exp, math
+
+func pow*(x: StUint, y: Natural): StUint {.inline.} =
+  when x.data is UintImpl:
+    result.data = x.data.pow(y)
+  else:
+    result.data = x.data ^ y

tests/property_based.nim

@@ -7,7 +7,7 @@
 #
 # at your option. This file may not be copied, modified, or distributed except according to those terms.
 
-import ../stint, unittest, quicktest
+import ../stint, unittest, quicktest, math
 
 const itercount = 1000
 
@@ -216,3 +216,16 @@ suite "Property-based testing (testing with random inputs) - uint64 on 64-bit /
         tz = tx div ty
 
     check(cast[uint](tz) == x div y)
+
+  quicktest "pow", itercount do(x: uint(min=0, max=hi), y: int(min = 0, max = high(int))):
+
+    when sizeof(int) == 8:
+      let
+        tx = cast[StUint[64]](x)
+        tz = tx.pow(y)
+    else:
+      let
+        tx = cast[StUint[32]](x)
+        tz = tx.pow(y)
+
+    check(cast[uint](tz) == x ^ y)

tests/property_based_uint256.nim

@@ -297,6 +297,8 @@ suite "Property-based testing (testing with random inputs) of Uint256":
       ttm_z = ttm_x mod ttm_y
       mp_z  = mp_x  mod mp_y
 
+    check ttm_z.asSt == mp_z
+
   quicktest "`div`", itercount do(x0: uint(min=0, max=hi),
                                 x1: uint(min=0, max=hi),
                                 x2: uint(min=0, max=hi),
@@ -319,3 +321,21 @@ suite "Property-based testing (testing with random inputs) of Uint256":
       ttm_z = ttm_x div ttm_y
       mp_z  = mp_x  div mp_y
 
+    check ttm_z.asSt == mp_z
+
+  quicktest "pow", itercount do(x0: uint(min=0, max=hi),
+                                x1: uint(min=0, max=hi),
+                                x2: uint(min=0, max=hi),
+                                x3: uint(min=0, max=hi),
+                                y : int(min=0, max=high(int))):
+
+    let
+      x = [x0, x1, x2, x3]
+      ttm_x = x.asTT
+      mp_x  = cast[StUint[256]](x)
+
+    let
+      ttm_z = ttm_x.pow(y.uint)
+      mp_z  = mp_x.pow y
+
+    check ttm_z.asSt == mp_z

tests/test_uint_exp.nim

@@ -0,0 +1,28 @@
+# Stint
+# Copyright 2018 Status Research & Development GmbH
+# Licensed under either of
+#
+#  * Apache License, version 2.0, ([LICENSE-APACHE](LICENSE-APACHE) or http://www.apache.org/licenses/LICENSE-2.0)
+#  * MIT license ([LICENSE-MIT](LICENSE-MIT) or http://opensource.org/licenses/MIT)
+#
+# at your option. This file may not be copied, modified, or distributed except according to those terms.
+
+import ../stint, unittest, math
+
+suite "Testing unsigned exponentiation":
+  test "Simple exponentiation 5^3":
+
+    let
+      a = 5'u64
+      b = 3
+      u = a.stuint(64)
+
+    check: cast[uint64](u.pow(b)) == a ^ b
+
+  test "12 ^ 34 == 4922235242952026704037113243122008064":
+    # https://www.wolframalpha.com/input/?i=12+%5E+34
+    let
+      a = 12.stuint(256)
+      b = 34
+
+    check: a.pow(b) == "4922235242952026704037113243122008064".u256