Thursday, January 21, 2010

Binary multiplication in Python

I came across a question on Stack OverFlow that references binary multiplication. This is a fun problem for a beginning Python coder.

Let's say we want to multiply 3 x 5 in binary:

  00101  ( 5)
x 00011 ( 3)
-------
00101 ( 5)
+ 01010 (10)
-------
01111 (15)


Let's talk generally by substituting m x n (so here m = 3 and n = 5).

No matter what the value of m, we will end up by adding together a series of bit-shifted versions of n, here: 00101 and 01010

5 = 00101
3 = 00011 = 2**0 + 2**1
3 x 5 = (5 << 0) + (5 << 1)


For any m, we require the set of exponents i, j, k.. such that:

2**i + 2**j + 2**k + ... = m

[This has to be done from the "top down." First find the largest power of 2 that is less than m, record that exponent, and take m = m % 2**e to continue... It's tricky, and that's why I took the easy way out, below.]

Now we sum

n << i
+ n << j
+ n << k
+ ... = m x n


It's very helpful that Python will do bit-shifting on ints:

>>> 5 << 1
10

I spent a while debugging code to generate the needed exponents, before I decided to cheat (a bit) by using binary representation:

>>> bin(3)
'0b11'

The first function below returns the exponents i,j,k...
The second does the bit-shifting. For efficiency, we do the factoring on the smaller of the two operands.

import random

def binary_exps(n):
L = list(bin(n)[2:])
L.reverse()
L = [i for i,c in enumerate(L) if c == '1']
return L

def multiply(x,y):
if y > x: x,y = y,x
result = 0
L = binary_exps(y)
for e in L:
result += x << e
return result, L

if __name__ == '__main__':
for i in range(1, 10):
n = 2**i-1
L = binary_exps(n)
print str(n).rjust(4), L
assert n == sum([2**e for e in L])
print

for i in range(20):
x = random.choice(range(1000))
y = random.choice(range(1000))
prod, L = multiply(x,y)
print str(x).rjust(4), str(y).rjust(4), L
assert prod == x * y

$ python multiply_binary.py 
1 [0]
3 [0, 1]
7 [0, 1, 2]
15 [0, 1, 2, 3]
31 [0, 1, 2, 3, 4]
63 [0, 1, 2, 3, 4, 5]
127 [0, 1, 2, 3, 4, 5, 6]
255 [0, 1, 2, 3, 4, 5, 6, 7]
511 [0, 1, 2, 3, 4, 5, 6, 7, 8]

406 819 [1, 2, 4, 7, 8]
902 817 [0, 4, 5, 8, 9]
755 115 [0, 1, 4, 5, 6]
963 170 [1, 3, 5, 7]
526 836 [1, 2, 3, 9]
236 228 [2, 5, 6, 7]
504 448 [6, 7, 8]
438 840 [1, 2, 4, 5, 7, 8]
225 791 [0, 5, 6, 7]
373 512 [0, 2, 4, 5, 6, 8]
307 230 [1, 2, 5, 6, 7]
672 565 [0, 2, 4, 5, 9]
389 137 [0, 3, 7]
444 406 [1, 2, 4, 7, 8]
587 707 [0, 1, 3, 6, 9]
412 246 [1, 2, 4, 5, 6, 7]
390 307 [0, 1, 4, 5, 8]
269 368 [0, 2, 3, 8]
146 553 [1, 4, 7]
90 182 [1, 3, 4, 6]