Python - 多桁乗算（Karatsuba 法）！

前回、 Python 3 で「標準（筆算）法」による多桁乗算アルゴリズムを実装する方法を紹介しました。

• Python - 多桁乗算（標準（筆算）法）！

今回は、「標準（筆算）法」より高速に乗算が可能な「Karatsuba 法」アルゴリズムを実装してみました。

0. 前提条件

• LMDE 2 (Linux Mint Debian Edition 2; 64bit) での作業を想定。
• Python 3.6.4 での作業を想定。
• 当方は他のバージョンとの共存環境であり、 python3.6, pip3.6 で 3.6 系を使用するようにしている。（適宜、置き換えて考えること）

1. アルゴリズムについて

当ブログ過去記事を参照。

• C++ - 多桁乗算（Karatsuba 法）！
• Ruby - 多桁乗算（Karatsuba 法）！

2. Python スクリプトの作成

• 敢えてオブジェクト指向で作成している。
• Shebang ストリング（1行目）では、フルパスでコマンド指定している。（当方の慣習）
• 数値計算ライブラリ NumPy を使用しない。（この程度の計算では、逆に2倍程度時間がかかってしまうため）
• 必要であれば、スクリプト内の定数を変更する。

File: multiply_karatsuba.py

#! /usr/local/bin/python3.6
"""
Multiplication of big-digit values with Karatsuba method
"""
import random
import sys
import traceback


class MultiplyKaratsuba:
    D_MAX = 1024  # Maximum number of digits (power of 2)
    D     = 1024  # Digits of computation (Multiples of 4 <= D_MAX)

    def __init__(self):
        self.a = [random.randrange(10) for _ in range(self.D)]
        self.b = [random.randrange(10) for _ in range(self.D)]

    def compute(self):
        """ Computation of multiplication """
        try:
            for i in range(self.D_MAX - len(self.a)):
                self.a.append(0)
            for i in range(self.D_MAX - len(self.b)):
                self.b.append(0)
            z = self.__multiply_karatsuba(self.a, self.b)
            z = self.__do_carry(z)
            self.__display(self.a, self.b, z)
        except Exception as e:
            raise

    def __multiply_normal(self, a, b):
        """ Normal multiplication

        :param  list a
        :param  list b
        :return list z
        """
        try:
            a_len, b_len = len(a), len(b)
            z = [0 for _ in range(a_len + b_len)]
            for j in range(b_len):
                for i in range(a_len):
                    z[j + i] += a[i] * b[j]
            return z
        except Exception as e:
            raise

    def __multiply_karatsuba(self, a, b):
        """ Karatsuba multiplication

        :param  list a
        :param  list b
        :return list z
        """
        try:
            t_len = len(a)
            if t_len <= 4:
                return self.__multiply_normal(a, b)
            a_0 = a[:(t_len // 2)]
            a_1 = a[(t_len // 2):]
            b_0 = b[:(t_len // 2)]
            b_1 = b[(t_len // 2):]
            v, w = [], []
            for i in range(t_len // 2):
                v.append(a_1[i] + a_0[i])
                w.append(b_1[i] + b_0[i])
            x_1 = self.__multiply_karatsuba(a_0, b_0)
            x_2 = self.__multiply_karatsuba(a_1, b_1)
            x_3 = self.__multiply_karatsuba(v,  w)
            for i in range(t_len):
                x_3[i] -= x_1[i] + x_2[i]
            z = x_1 + x_2
            for i in range(t_len):
                z[i + t_len//2] += x_3[i]
            return z
        except Exception as e:
            raise

    def __do_carry(self, a):
        """ Process of carrying

        :param  list a
        :return list a
        """
        cr = 0
        try:
            for i in range(len(a)):
                a[i] += cr
                cr = a[i] // 10
                a[i] -= cr * 10
            if cr != 0:
                print("[ OVERFLOW!! ] ", cr)
            return a
        except Exception as e:
            raise

    def __display(self, a, b, z):
        """ Display

        :param list a
        :param list b
        :param list z
        """
        a_len = self.D_MAX
        b_len = self.D_MAX
        z_len = self.D_MAX * 2
        try:
            while a[a_len - 1] == 0:
                if a[a_len - 1] == 0:
                    a_len -= 1
            while b[b_len - 1] == 0:
                if b[b_len - 1] == 0:
                    b_len -= 1
            while z[z_len - 1] == 0:
                if z[z_len - 1] == 0:
                    z_len -= 1
            print("a =")
            for i in reversed(range(a_len)):
                print(a[i], end="")
                if (a_len - i) % 10 == 0:
                    print(" ", end="")
                if (a_len - i) % 50 == 0:
                    print()
            print()
            print("b =")
            for i in reversed(range(b_len)):
                print(b[i], end="")
                if (b_len - i) % 10 == 0:
                    print(" ", end="")
                if (b_len - i) % 50 == 0:
                    print()
            print()
            print("z =")
            for i in reversed(range(z_len)):
                print(z[i], end="")
                if (z_len - i) % 10 == 0:
                    print(" ", end="")
                if (z_len - i) % 50 == 0:
                    print()
            print()
        except Exception as e:
            raise


if __name__ == '__main__':
    try:
        obj = MultiplyKaratsuba()
        obj.compute()
    except Exception as e:
        traceback.print_exc()
        sys.exit(1)

• Gist - Python script to multiply big-digit values with Karatsuba method.

3. Python スクリプトの実行

まず、実行権限を付与。

$ chmod +x multiply_karatsuba.py

そして、実行。

$ ./multiply_karatsuba.py
a =
2105985447 9026926375 6163318067 1521930549 6102825804
1732057778 1611724768 2908201805 0244698557 0709918605
3288302496 4685317623 4898190784 4604144986 3284474952
3815865438 2836148231 2276177859 5371315395 9309499404
0566467922 5640804557 4467398356 7424949681 5178552540
0541446090 9708451086 7487095916 4646885025 2741712589
2363632763 0729423690 1494709138 6427545614 3028890118
4889790537 0353275730 9503292442 5163455511 5505731441
8592814590 9017490715 7015863874 1973673448 7936782139
8928318070 9515820526 3695226727 0786268578 8720379890
6887133727 2018526408 6131597125 6321843839 0807036026
3769160191 5458712924 4834971170 0441452910 5739231133
0670340821 6600754190 0522297563 8211057443 1231388408
3316296015 0459410484 8764822388 7751806458 6495800912
5564499453 2656329959 3722669651 7720216391 6356637691
3337369251 3648011916 9573403651 3886387186 2826292213
7560755759 9805658710 3767623191 4379788385 7168928193
3975006892 7048951344 3313637262 4872387200 2683660287
0676925624 9068704392 8887230169 6092366260 3809264412
8164893414 7875628093 7481381683 0046583990 3655373497
1748152522 2338977371 6948
b =
4366469670 3352625316 5510628732 7279216713 9157308666
8404530809 5155246561 9742812159 3712625456 1332100320
1707557823 4106256438 0173900781 3010030878 4575894860
4367607055 8736235625 9100602978 1673553804 2210299395
6686456288 5410955347 2188689413 4870874031 9921652732
9578669947 6155568279 4319874491 0773822156 4283762043
3744141415 6586469207 8334719319 9332699957 1296982503
3013589754 4141518291 8395808984 7573901094 8772990979
3589667620 4966730949 4557350864 1346004040 1666278583
8383072207 6458936430 9703269027 9983478606 6353194116
6839882498 7202316199 3088389089 2028400979 5573595533
0872930081 0735540881 0101544157 6600790471 2281515618
4515107231 9208853815 0984860225 5179140820 9457673331
2100093633 8850378825 1781382974 4695469940 5116060013
4142133163 6491857677 6398988539 5696980766 3479164616
4820409853 3521017829 1941410440 0677222815 9608160501
4125204122 9518267483 5078059055 0055299138 9173535281
0211731935 7691628736 5328243792 6164518594 8050934433
7320800390 5109094720 4205090660 0997232578 3995922913
0436793066 6125787094 3388887009 7226333094 0070116097
3124359899 1053911781 8075
z =
9195721584 4345305261 3784372239 2619730025 9709288555
0625891700 0425091970 5413580730 1838410598 1164579193
3556857349 6319147194 0405916033 5584528780 3340215384
6800763238 4643269282 1879868418 2391971025 8849100804
6633941047 1304637342 9062754152 2987910268 5053507656
7936530947 5464952839 5584313960 4493822379 2041595745
3207454440 3097602371 8089400555 9605381499 2253277639
8266368350 2608566970 4464472909 7999555769 2708536577
0616719377 0555861205 9252613165 8838805388 6785051192
6607613096 9015161138 8764149133 5742091218 3785052342
4937468039 0912853664 3108498441 5448459921 0240328721
4716227000 3334750768 5660325424 2343585147 2633224213
2134953985 9067375586 9519732307 7708998788 5112124164
8732422570 6197477692 9538435683 8823915544 9517333936
9274932251 9721121276 0302054521 7249457306 8634514334
2964964063 4858530853 0401493135 1282742659 0080652728
5397477551 0094744816 6624970112 6483824369 4343913667
8386831589 8791180640 4444539898 5200149392 7861376145
0401281911 8056943064 6130910692 9125515495 3741340417
7047486957 9477989920 0474131121 1210178329 4417738382
1269471393 8180898412 9881742528 6382100174 0883434176
3236106165 9341657226 1629754331 1332530711 5862296183
9501218376 8186735690 9846766197 3758818963 5915329550
4221773538 8288527196 9462799411 7394177062 6672912386
6336959888 1816617528 7817188649 3153015340 7828861162
3493890190 7564529087 3756497382 7667478003 8686526841
5611259682 9783562867 1923519039 1203474353 2247999689
0735672946 9842721936 3639217701 9643387694 9284007798
2655579913 7146465504 5695644122 1017911119 1004716200
9478455686 2663063591 6358475221 0014850178 7072237219
2505209139 6934611505 2539130584 3486847403 8666463342
7267992432 9059939037 5244690943 0276464594 8658846549
3550855848 3233476890 4590488286 3377284029 8654472183
1156295479 8808632082 5449318945 7152398719 4788700237
9456043047 5812881514 5675418478 0060413341 8269359936
3611305317 3675808350 5630467111 5673108244 3215605764
1592213906 3322057143 9017160501 5968667798 8667476788
4363142935 2835338163 4639119738 7725710450 5582346635
0664304068 3422461938 3090967032 1825513189 3488848468
2542353714 0881370058 1739983017 4332759708 1969165524
5022481864 0192757494 6014532112 8122898840 8235100

---

以上