@@ -70,3 +70,11 @@ euclide a b = let (d', u', v') = euclide b (mod a b) in (d', v', u' - (div a b)
modInv e n = let (_, d, _) = euclide e n in d
+-- Return x ^ k (mod) n
+expMod x k n =
+ if k == 0
+ then 1
+ else if even k
+ then expMod (mod ((mod x n) * (mod x n)) n) (div k 2) n
+ else (mod x n) * expMod ((mod x n) * (mod x n)) (div k 2) n
+