整數的故事(4)——Karastuba算法

　　我們在小學就學過用豎式計算兩個多位數的乘法：

　　這個過程簡單而繁瑣，沒有最強大腦的普通大眾通常是用計算器代替的。然而對於超大整數的乘法，計算器也未必靠得住，它還存在“溢出”一說。這就需要我們自行編寫算法了。

豎式算法

　　雖然對於Python來說，不必太過關心整數的長度和溢出問題，但對於其它編程語言就未必了。這里我們暫且拋開語言本身的特性，只關注算法本身。假設輸入的兩個長整數x和y，它們的乘積將會溢出，所以需要將乘積轉換成字符串，根據乘法豎式的運算規則，很容易寫出下面的代碼：

# 豎式法計算兩個整數相乘
def multi(x, y):
    x_str = str(x)
    y_str = str(y)
    x_len = len(x_str)
    y_len = len(y_str)
    z = [0] * (x_len + y_len)

    for i in range(x_len):
        for j in range(y_len):
            z[x_len - i - 1 + y_len - j - 1] += int(x_str[i]) * int(y_str[j])

    for i in range(len(z) - 1):
        if z[i] >= 10:
            z[i + 1] += (int(z[i]) // 10)
            z[i] = int(z[i]) % 10

    return list2str(z)

# 將z轉換成字符串並刪除左側的0
def list2str(z):
    result = [str(i) for i in z]
    return ''.join(result[::-1]).lstrip('0')

# 打印運行結果
def paint(a, b):
    print('{0} * {1} = {2}'.format(a, b, multi(a, b)))

if __name__ == '__main__':
    paint(123,321)
    paint(123,456)
    paint(123456789000, 987654321000)

　　代碼中9~11行用兩個循環模擬了乘法計算的過程，以123×321為例，在循環結束后z將存儲下面的數據：

　　11行的for循環是處理進位問題。最后將列表轉換為字符串，再去掉多余的0，打印結果：

Karastuba算法

　　豎式乘法偏向於使用蠻力，Karastuba博士在1960年提出了一個更簡單的算法，其思想是把兩個大整數的乘法轉化為若干次小規模的乘法和少量的加法，這就是Karastuba算法。

　　對於兩個n位的大整數x和y，可以把x和y分解成兩部分：

　　例如：

　　是不是有點似成相識？沒錯，這實際上是利用了歐幾里德算式將一個整數分解成m=qn+r的形式。現在x和y的乘積可以表示為：

　　這就把原來的大整數乘法變成了四次效較小規模的乘法（其中10ⁿ的運算可以通過位移高效處理）和少量加法。上式還可以更進一步：

　　看起來更復雜了，但是對於計算機來說，x₁y₁和x₀y₀已經計算過了，不需要再次計算。x₁y₀+x₀y₁被轉換成一次乘法和少量的加法，多一個加法運算對時間復雜度沒有影響，而減少一個乘法卻能減少時間復雜度。對每一個乘法都進行類似的分解，反復迭代x_iy_i，直到其中一個乘數只有1位為止。按照這種思路可以編寫新的乘法運算代碼：

# karastuba算法計算兩個n位的大整數乘法, x >=0, y >= 0
def karastuba(x, y, n):
    if x == 0 or y == 0:
        return 0
    elif n == 1:
        return x * y

    k = n // 2
    x1 = x // (10 ** k)
    x0 = x % (10 ** k)
    y1 = y // (10 ** k)
    y0 = y % (10 ** k)
    z0 = karastuba(x0, y0, k)         # 計算x0y0
    z1 = karastuba(x1, y1, k)         # 計算x1y1
    z2 = karastuba((x1 + x0), (y0 + y1), k) - z1 - z0

    return z1 * (10 ** n) + z2 * (10 ** k) + z0

　　然而運行時會發現這段代碼很難生效，原因是計算時要求的環境太過理想——每次迭代時x_i和y_i的位數都必須相同。這就需要重新審視Karastuba算法，看看非理想狀態下是如何計算的。

　　假設x和y分別是m位和n位的大整數，x和y可以這樣分解：

　　反復迭代x_iy_i，直到其中一個乘數只有1位為止。

　　示例： 123×321 = ?

　　

 

　　現在可以編寫能夠正確運行的大整數乘法代碼：

def karastuba(x, y):
    ''' karastuba算法計算兩個n位的大整數乘法, x >=0, y >= 0 '''
    if x == 0 or y == 0:
        return 0
    m, n = len(str(x)), len(str(y)) # x和y的位數
    if m == 1 or n == 1: # 如果x或y只有1位，直接計算結果
        return x * y
    m //= 2
    x1, x0 = x // (10 ** m), x % (10 ** m) # 分解x
    n //= 2
    y1, y0 = y // (10 ** n), y % (10 ** n) # 分解y
    # 迭代分解夠的4個較小規模的乘法
    x1y1 = karastuba(x1, y1)
    x1y0 = karastuba(x1, y0)
    x0y1 = karastuba(x0, y1)
    x0y0 = karastuba(x0, y0)
    return x1y1 * (10 ** (m + n)) + x1y0 * (10 ** m) + x0y1 * (10 ** n) + x0y0

def paint(x, y):
    ''' 在控制台打印karastuba(x, y)的運行結果 '''
    print('{0} * {1} = {2}, ({3})'.format(x, y, karastuba(x, y), x * y))

if __name__ == '__main__':
    paint(123, 321)
    paint(123456789, 987456)
    paint(1234567891234567, 1234567891234567)

 

　　先看看windows計算器下1234567891234567×1234567891234567的運行結果：

　　計算器已經無法給出精確的結果，但karastuba沒有問題：

　　對於非10進制整數，Karastuba算法依然適用。

　　

---

　　 作者：我是8位的

　　出處：http://www.cnblogs.com/bigmonkey

　　本文以學習、研究和分享為主，如需轉載，請聯系本人，標明作者和出處，非商業用途！ 

　　掃描二維碼關注公眾號“我是8位的”