In [1]:
import numpy as np # used for generating random numbers
In [2]:
def int_to_big(x):
if x == 0:
return [0]
z = []
while x > 0:
t = x % 10
z.append(t)
x //= 10
trim(z)
return z
def big_to_int(x):
z, p = 0, 1
for d in x:
z += p * d
p *= 10
return z
In [3]:
from itertools import zip_longest
def trim(z):
while len(z) > 1 and z[-1] == 0:
z.pop(-1)
def add(x, y):
z, carry = [], 0
for r, s in zip_longest(x, y, fillvalue=0):
carry += r + s
z.append(carry % 10)
carry //= 10
if carry:
z.append(carry)
return z
def subtract(x, y):
z, carry = [], 0
for r, s in zip_longest(x, y, fillvalue=0):
carry += r - s
z.append(carry % 10)
carry //= 10
trim(z)
return z
In [4]:
def karatsuba(x, y):
# ensure same length
while len(x) < len(y):
x.append(0)
while len(x) > len(y):
y.append(0)
# length
n = len(x)
half = n // 2
if n == 1:
return add([x[0] * y[0]], [])
# cut-off for improved efficiency
if n <= 50:
a = big_to_int(x)
b = big_to_int(y)
z = a * b
return int_to_big(z)
x0, x1 = x[:half], x[half:]
y0, y1 = y[:half], y[half:]
# x = x0x1
# y = y0y1
# z0 = x0 * y0
# z1 = x1 * y1
# z2 = (x0 + x1) * (y0 + y1)
# z2 = z2 - (z0 + z1)
z0 = karatsuba(x0, y0)
z1 = karatsuba(x1, y1)
z2 = karatsuba(add(x0, x1), add(y0, y1))
z2 = subtract(z2, add(z0, z1))
z = add(z0, [0] * (half << 1) + z1)
z = add(z, [0] * half + z2)
return z
In [5]:
def multiply(x, y):
xb = int_to_big(x)
yb = int_to_big(y)
zb = karatsuba(xb, yb)
return big_to_int(zb)
def test(x, y):
z = multiply(x, y)
assert x * y == z
print("{} x {} = {}".format(x, y, z))
In [6]:
def gen_long(n):
x = ''.join(map(str, np.random.randint(0, 10, n)))
return int(x)
In [7]:
test(1432423423420, 12321312332131233)
test(8931283129323420, 1233123602345430533)
In [8]:
tests = 30
for _ in range(tests):
n = np.random.randint(1, 15)
x, y = gen_long(n), gen_long(n)
test(int(x), int(y))
In [9]:
%%time
a, b = gen_long(1000), gen_long(1000)
z = multiply(a, b)
assert z == a * b
In [10]:
%%time
a, b = gen_long(20000), gen_long(20000)
z = multiply(a, b)
assert z == a * b
Karatsuba's algorithm is already implemented in Python. Check this package.
In [11]:
from karatsuba import *
In [12]:
def power_of_two(x):
p = 1
while p < x:
p <<= 1
return p
def reverse(num):
return int(str(num)[::-1])
In [13]:
def kat_multiply(x, y):
if x == 0 or y == 0:
return 0
xs = list(map(int, str(x)))
ys = list(map(int, str(y)))
n = power_of_two(max(len(xs), len(ys)))
plan = make_plan(range(n), range(n))
xs = [0] * (n - len(xs)) + xs
ys = [0] * (n - len(ys)) + ys
zs = plan(xs, ys)
zs.pop(-1)
zs = list(reversed(zs))
while zs[-1] == 0:
zs.pop(-1)
ans = 0
for p, d in enumerate(zs):
ans += d * 10 ** p
return ans
In [14]:
tests = 30
for _ in range(tests):
n = np.random.randint(1, 15)
x, y = gen_long(n), gen_long(n)
z = kat_multiply(x, y)
assert z == x * y
print("{} x {} = {}".format(x, y, z))
In [15]:
%%time
a, b = gen_long(100), gen_long(100)
z = kat_multiply(a, b)
assert z == a * b