pythoniku
pythoniku

Reputation: 3662

Karatsuba algorithm incorrect result

I just simply followed the pseudo code on wiki http://en.wikipedia.org/wiki/Karatsuba_algorithm But the result of this implementation is very unstable. It works sometimes but in case like 100*100. It does fail. What I missed here? please take a look.

from math import *
f = lambda x: (int(x) & 1 and True) and 1
def fast_multiply( x = "100", y = "100"):
    print "input "+x+" | "+y
    int_buff = map( int, [x, y])
    if int_buff[0] < 10 or int_buff[1] < 10:
        #print "lol"
        return int_buff[0]*int_buff[1]

    degree = max( x.__len__(), y.__len__())

    higher_x, lower_x = x[ : int( ceil( len(x) / 2.0))], x[ len(x)/2 +f(len(x)):]
    higher_y, lower_y = y[ : int( ceil( len(y) / 2.0))], y[ len(y)/2 +f(len(y)):]
    #print lower_x+" & "+lower_y
    z0 = fast_multiply(lower_x, lower_y) #z0 = 0
    z1 = fast_multiply(str(int(lower_x)+int(higher_x)), str(int(lower_y)+int(higher_y)))
    z2 = fast_multiply(higher_x, higher_y)
    print "debug "+str(z0)+" "+str(z1)+" "+str(z2)
    return z2*(10**degree) + (z1-z2-z0)*(10**(degree/2))+z0




if __name__ == '__main__':
    print fast_multiply()

I have noticed in the case 100*100 z2 will be 100 which is correct. This gives z2*(10**3)=100000 which is definitely wrong...

Upvotes: 0

Views: 788

Answers (1)

Richard
Richard

Reputation: 362

The pseudocode you used was wrong. The problem was in z2*(10**degree). You should have raised the base to 2*m where m is what you meant to calculate with int( ceil(len(x) / 2.0)) (len(x) and len(y) should both have been degree).

I couldn't resist refactoring it... a little. I used the names from the definitions on the wiki. It would be straightforward to implement it with an arbitrary base, but I stuck with 10 for simplicity.

def kmult(x, y):
    if min(x, y) < 10:
        return x * y

    m = half_ceil(degree(max(x, y)))

    x1, x0 = decompose(x, m)
    y1, y0 = decompose(y, m)

    z2 = kmult(x1, y1)
    z0 = kmult(x0, y0)
    z1 = kmult(x1 + x0, y1 + y0) - z2 - z0

    xy = z2 * 10**(2*m)  +  z1 * 10**m  +  z0
    return xy


def decompose(x, m):
    return x // 10 ** m, x % 10 ** m

def degree(x):
    return len(str(x))

def half_ceil(n):
    return n // 2 + (n & 1)

Testing:

print kmult(100, 100)

def test_kmult(r):
    for x, y in [(a, b) for b in range(r+1) for a in range(r+1)]:
        if kmult(x, y) != x * y:
            print('fail')
            break
    else:
        print('success')


test_kmult(100)

Result:

10000
success

Upvotes: 2

Related Questions