乘法算法

乘法算法是計算兩個數值相乘乘積的演算法。為了提高運算效率，不同大小的數字適用不同的乘法算法。自十進位數字系統誕生以來，就已開始發展乘法算法。

## 網格法

網格式乘法 （或盒式法）是一種用於給小學生進行乘法計算啟蒙的簡單乘法算法。自1990年代後期以來，它一直是英格蘭和威爾斯公立小學數學課程的標準部分。

將兩個乘數分解（「劃分」）為百位、十位或個位的部分，然後在相對簡單的乘法步驟中直接算出這些部分的乘積，再將其一一相加以算出最終結果。用這種方法計算 $34 \times 13$，可以畫出如下網格：

```
  300
   40
   90
 + 12
 ————
  442
```

×｜30｜4

10｜300｜40

3｜90｜12

最後再通過一次求和或逐行累加得到結果（見右），實質是計算了 $(300 + 40) + (90 + 12) = 340 + 102 = 442$。

這種計算方法（儘管不一定顯式地用網格排列乘數）也稱為「部分乘積演算法」。它的本質是分別計算簡單的個位乘法，所有加法都留到最後的「加總」階段。

儘管隨著位數的增加，伴隨的中間過程越來越繁瑣，但網格方法原則上可以應用於任何大小的乘數。一般認為此方式是引入多位乘法概念時，有用的顯式方法，而在當今利用計算器或電子表格就能完成大多數計算的時代，它可能是某些學生唯一需要的乘法算法。

## 長乘法

位值制數字系統在人類社會的壟斷性地位，使得長乘法成為最自然的乘法算法，世界各地的學校都會教學生使用這種演算法完成日常的乘法計算。其內容是：將一個乘數乘以另一個乘數的每位數字，然後將所有數字按照適當的位置相加得出結果。這需要預先記住所有個位數字兩兩相乘後的結果，即華人世界常見的乘法表。長乘法又被稱為教科書乘法或標準乘法。

若在十進制下，人工計算較大的數字加乘，長乘法是常見的演算法。電腦在二進制下，也有類似的演算法「移位和相加」，不少現代的處理器已設計最佳化的乘法線路，以進行快速的乘法計算，但其代價是硬體會變的更加複雜。若是在紙上計算長乘法，會先將所有數字寫下，然後相加，若是使用珠算，會在計算完每一個乘積後，直接和已計算的和加總。

### 示例

此例用長乘法計算23,958,233（被乘數）乘以5,830（乘數），結果是139,676,498,390（積）。

```
        23958233
  ×         5830
  ———————————————
        00000000 ( =      23,958,233 ×     0)
       71874699  ( =      23,958,233 ×    30)
     191665864   ( =      23,958,233 ×   800)
  + 119791165    ( =      23,958,233 × 5,000)
  ———————————————
    139676498390 ( = 139,676,498,390        )
```

以下的虛擬碼描述上述乘法的步驟。此作法會持續的更新乘積，在乘完所有數字後，乘積即為結果。其中+=運算子用來表示將某變數原有的值加上其他的值，再存回該變數中（類似在Java及C語言中的用法），以讓程式緊湊。

```pascal
multiply(a[1..p], b[1..q], base)                            // index 1對應數字的個位數
  product = [1..p+q]                                        // 先預留乘積的空間
  for b_i = 1 to q                                          // 針對b的每一位數
    carry = 0
    for a_i = 1 to p                                        // 針對a的每一位數
      product[a_i + b_i - 1] += carry + a[a_i] * b[b_i]
      carry = product[a_i + b_i - 1] / base
      product[a_i + b_i - 1] = product[a_i + b_i - 1] mod base
    product[b_i + p] = carry                               // 最高位數是最後一次計算乘法的進位
  return product
```

### 優化空間複雜度

假設在基數為$D$，令 $n$ 為兩個輸入數字位數長度的總和，若其結果需保存在存儲器中，則其空間複雜度通常為$\Theta(n)$。然而在某些應用中，不需要把整個結果保存在記憶體中，而可以在計算時將數字以字元流的方式輸出（例如輸出到系統終端或文件中）。在這些情況下，長乘法的優勢在於可以將其輕鬆表示為對數空間演算法——也就是說，只需要一個工作空間與輸入中位數的對數成正比的演算法（$\Theta(\log n)$），這是乘數之積的雙重對數（$\log \log N$）。注意，乘數本身仍保留在記憶體中，並且在此分析中忽略其 $\Theta(n)$ 的空間。

只要知道之前乘積產生的進位，就可以由最低位起，一位一位的計算乘積的每一位數字，這是優化空間複雜度的關鍵。令ai和bi分別是被乘數及乘數的第i位數字，a和 b的最高位都已補0，補到n位數，而ri是乘積的第i位，而ci是計算ri後產生的進位（i=1表示是個位），則

$$\begin{align} r_i &= \left( c_{i-1} + \sum_{j+k=i+1} a_j b_k \right) \mod D \\ c_i &= \left\lfloor (c_{i-1} + \sum_{j+k=i} a_j b_k) / D \right\rfloor \\ c_0 &= 0 \end{align}$$

or

$$c_i = \left( \sum_{m = 0}^{i - 2} \sum_{j + k = i - m} a_j b_k \right) / D.$$

- 簡單的推論可以證明進位不會超過n，ri的總和不會超過D
- n：第一欄的進位是0，其他欄最多有n個數字，最多也只會從前一欄進位n。總和最多是D
- n，進位到下一位的最多是D
- n / D或n。因此這些數字都可以用O(log n)位數儲存。

偽代碼下的對數空間演演算法為

```pascal
multiply(a[1..p], b[1..q], base)                  // Operands containing rightmost digits at index 1
    tot = 0
    for ri = 1 to p + q - 1                       //For each digit of result
        for bi = MAX(1, ri - p + 1) to MIN(ri, q) //Digits from b that need to be considered
            ai = ri − bi + 1                      //Digits from a follow "symmetry"
            tot = tot + (a[ai] * b[bi])
        product[ri] = tot mod base
        tot = floor(tot / base)
    product[p+q] = tot mod base                   //Last digit of the result comes from last carry
    return product
```

### 在電腦中的用法

有些積體電路會實現乘法，可能是用硬體或是微程序，可能針對各種整數及浮點數。在高精度計算中，在乘比較小的數字時，用底數為2w的長乘法，其中w為字元的位元數。

若要用此方式乘二個n位數數字，約需要n2次的運算。若用比較型式的方式表示，用數字位數這個自然的數字大小度量，用長乘法乘二個n位數數字的時間複雜度是Θ(n2)。

若要用軟體實現，長乘法演演算法需要處理在加法時會出現的溢位問題，可能會很耗成本。典型的作法是用比較小的基底來進行運算，例如b，而讓8b為機器中表示的整數大小。這樣可以進行幾次運算之後才會溢位。若數字太大，可以只加數字的一部份到結果中，或是進位，將剩下的數字映射為一個小於b的較小數字。此作法稱為「正規化」。Richard Brent有在其Fortran軟體包MP中使用此一作法。

## 格子乘法

<image>

<image>

格子乘法在演演算法上和長乘法等效。需要在紙上劃格子，以引導計算，並且將乘法和加法分為二個步驟進行。這是在1202年在斐波那契的《計算書》中引進到歐洲。斐波那契用此方式心算，用左手和右手處理計算的中間值。16世紀的馬特拉克·那西在《Umdet-ul Hisab》書籍上列出了此演演算法六種不同的變體。在奧斯曼帝國的恩德倫學校中，廣為使用此一演演算法。納皮爾的骨頭（Napier's bones）或稱為納皮爾算籌（Napier's rods）也是用此方式，是約翰·納皮爾過世的1617年發表。

如圖所示，被乘數和乘數寫在格子的上方和右邊，此方法在花拉子米的《算數》（Arithmetic）書中有提到，這是斐波那契的來源之一，曾被Sigler在2002年的《Fibonacci's Liber Abaci》中提到。

- 在乘法階段，格子中會填上其對應行和列上所寫，被乘數和乘數的位數相乘後的二位數乘積：格子中左上方的三角形會列十位數字，格子中右下方的三角形會列個位數字。
- 在加法階段，會延著斜線將數字相加，將相加後的個位數寫在格子的最下方，對應斜線，若有進位，將進位寫在下一個斜線區域格子的右方或上方。
- 下方的數字即為兩數字相乘的積。

### 示例

考慮以下較複雜的計算，考慮23,958,233乘以5,830，結果是139,676,498,390。不過以下計算中，其中的23,958,233是寫在格子的上方，5,830寫在右方。將乘積寫在格子的二個三角形中，將和寫在左邊和下方。結果條列如下：

```
     2   3   9   5   8   2   3   3
   +---+---+---+---+---+---+---+---+-
   |1 /|1 /|4 /|2 /|4 /|1 /|1 /|1 /|
   | / | / | / | / | / | / | / | / | 5
 01|/ 0|/ 5|/ 5|/ 5|/ 0|/ 0|/ 5|/ 5|
   +---+---+---+---+---+---+---+---+-
   |1 /|2 /|7 /|4 /|6 /|1 /|2 /|2 /|
   | / | / | / | / | / | / | / | / | 8
 02|/ 6|/ 4|/ 2|/ 0|/ 4|/ 6|/ 4|/ 4|
   +---+---+---+---+---+---+---+---+-
   |0 /|0 /|2 /|1 /|2 /|0 /|0 /|0 /|
   | / | / | / | / | / | / | / | / | 3
 17|/ 6|/ 9|/ 7|/ 5|/ 4|/ 6|/ 9|/ 9|
   +---+---+---+---+---+---+---+---+-
   |0 /|0 /|0 /|0 /|0 /|0 /|0 /|0 /|
   | / | / | / | / | / | / | / | / | 0
 24|/ 0|/ 0|/ 0|/ 0|/ 0|/ 0|/ 0|/ 0|
   +---+---+---+---+---+---+---+---+-
     26  15  13  18  17  13  09  00
```

```
 01
 002
 0017
 00024
 000026
 0000015
 00000013
 000000018
 0000000017
 00000000013
 000000000009
 0000000000000
 —————————————
  139676498390
```

= 139,676,498,390

## 二進制或農夫演演算法

農夫演演算法是廣為在農民中使用的演算法，其中不需要像長乘法一樣記憶乘法表。此演演算法曾在古埃及使用，農夫演演算法的優點是可以快速教授、不需要記憶、若沒有紙筆的話，也可以用其他東西輔助計算（例如籌碼），缺點是步驟比長除法要長，在計算大數字的乘法時不方便。

### 說明

在紙上，寫下被乘數，被乘數的下方寫下被乘數反覆除二（且無條件捨去小數）的結果，被乘數的旁邊寫下乘數，乘數下方寫上乘數反覆乘二的結果。一直到被乘數那一欄為1為止。將乘數那一欄的數字相加，但若被乘數那一欄為偶數，跳過此數字不累加，所得結果即為乘積。

### 示例

以下用農夫演演算法計算11乘以3。

```
十進制:     二進制:
11   3       1011  11
5    6       101  110
2   
12
       10  
1100

1   24       1  11000
    ——         ——————
    33         100001
```

說明上述的步驟：

- 二欄的最上方分別是11和3。
- 11除以2（5.5）無條件捨去，變成5，3乘以2（6）。
- 5除以2（2.5）無條件捨去，變成2，6乘以2（12），這一行的最左邊（2）是偶數，因此12要劃掉，不能累加到乘積中。
- 2除以2（1），12乘以2（24）
- 將沒有劃掉的數字相加：3 + 6 + 24 = 33.

此作法有效的原因是因為乘法滿足分配律，因此：

$$\begin{align} 3 \times 11 & = 3 \times (1\times 2^0 + 1\times 2^1 + 0\times 2^2 + 1\times 2^3) \\ & = 3 \times (1 + 2 + 8) \\ & = 3 + 6 + 24 \\ & = 33. \end{align}$$

以下是一個複雜的例子，仍用之前用過的範例（23,958,233 and 5,830）：

```
Decimal:             Binary:
5830  
23958233
       1011011000110  
1011011011001001011011001

2915  47916466       101101100011  10110110110010010110110010
1457  95832932       10110110001  101101101100100101101100100
728  
191665864
       1011011000  
1011011011001001011011001000

364  
383331728
       101101100  
10110110110010010110110010000

182  
766663456
       10110110  
101101101100100101101100100000

91  1533326912       1011011  1011011011001001011011001000000
45  3066653824       101101  10110110110010010110110010000000
22  
6133307648
       10110  
101101101100100101101100100000000

11 12266615296       1011  1011011011001001011011001000000000
5  24533230592       101  10110110110010010110110010000000000
2  
49066461184
       10  
101101101100100101101100100000000000

1  98132922368       1  <u>1011011011001001011011001000000000000</u>
  ————————————          1022143253354344244353353243222210110 （進位前）
  139676498390         10000010000101010111100011100111010110
```

## 電腦中的二進制乘法

這是農夫演演算法的變體。

二進制下，長乘法會變成很簡單的處理，針對每一個被乘數位數的1，將乘數往左位移適當的位元數，將乘數各次位移的結果相加，即為乘積。在一些中央處理器中，加法和位移的速度會比乘法要快，尤其被乘數不大，或是幾乎是定值的情形。

## 移位和相加

以往的電腦會用「移位和相加」的乘法來計算小整數的乘法。二進制的長乘法和農夫演演算法都可以簡化到相同的演演算法。

在二進制下，將數字和二進制的一個位元相乘，可以簡化為各位元的邏輯與運算。會將算出來的部份和相加。大部份目前的微處理器都會以此方式或是其他類似方式（例如布斯乘法算法）實現不同整數以及浮點長度的乘法，可能是用乘法器中或是微程序。

目前的處理器，位元的移位運算會比乘法要快很多，因此可以用往左移位代替乘以2的次冪。乘以常數的乘法也可以用一連串的移位和加減法來代替。例如，以下是只有移位及相加來表示乘以10的演演算法。

```
((x << 2) + x) << 1 # 此處是用 (x*2^2 + x)*2 來計算 10*x
(x << 3) + (x << 1) # 此處是用 x*2^3 + x*2 來計算 10*x
```

有時的移位和加減法的組合會比硬體乘法器還快。

若電腦沒有乘法的運算，需要用上述方式來進行計算，若將向左移位表示為相同的數加二次（兩者在邏輯上等價），也可以延伸到浮點數。

## 平方的四分之一乘法

可以用平方的算式，計算二個數的乘積，有些來源的算式中會包括取整函數，此方式源自於巴比倫數學（西元前2000-1600年）。

$$\left\lfloor \frac{\left(x+y\right)^2}{4} \right\rfloor - \left\lfloor \frac{\left(x-y\right)^2}{4} \right\rfloor = \frac{1}{4}\left(\left(x^2+2xy+y^2\right) - \left(x^2-2xy+y^2\right)\right) = \frac{1}{4}\left(4xy\right) = xy.$$

x和y為整數，若x+y和x−y中有一個為奇數，另一個也一定會是奇數，因此上式中若有一項在取整數之前有分數，另一項也會有，兩者會相消，不會影響所得的結果。以下是從0到18計算平方的四分之一取整數的對照表，可以計算9×9以內的乘法：

n ｜ 0｜ 1｜ 2｜ 3｜ 4｜ 5｜ 6｜7｜8｜9｜10｜11｜12｜13｜14｜15｜16｜17｜18

⌊n2/4⌋｜0｜0｜1｜2｜4｜6｜9｜12｜16｜20｜25｜30｜36｜42｜49｜56｜64｜72｜81

若要計算9乘3，其和以及差分別是12以及6，根據上式可得36和9，兩者的差是27，這就是9乘3的乘積。

Antoine Voisin在1817年有出版一個平方的四分之一的表，從1列到1000，可以幫助乘法的運算。Samuel Laundy在1856年有出版過1到100000的表，Joseph Blater在1888年出版了1到200000的表。

平方的四分之一乘法曾用在模擬計算機上，以形成二訊號相乘的模擬信號。二個輸入電壓的和或差可以用運算放大器達到。電壓平方可以用分段線性函數電路達成。最後二個平方的四分之一以及兩電壓的差以及再用運算放大器實現。

Everett L. Johnson在1980年時提出將此方式應用在數碼乘法器中。為了產生8位元整數的乘積，數位裝置計算兩數的和以及差，查表計算平方，計算差值，再往右移位二個位元。

平方的四分之一乘法對8位元的系統有益，可以在沒有硬體乘法器的情形下實現乘法。Charles Putney有將此演演算法實現在MOS 6502中。

## 複數乘法演演算法

複數乘法包括四個實數乘法以及二個加法。

$$(a+bi) (c+di) = (ac-bd) + (bc+ad)i.$$

或

$$\begin{array}{c|c|c} \times & a & bi \\ \hline c & ac & bci \\ \hline di & adi & -bd \end{array}$$

不過有辦法將實數乘法由四個減少到三個

乘積(a + bi) · (c + di)可以用以下方式計算。

k1 = c · (a + b)

k2 = a · (d − c)

k3 = b · (c + d)

實部 = k1 − k3

虛部 = k1 + k2.

此演演算法只用了三個乘法，比原來的四個乘法少一個，但加減法由原來的二個增加到五個。假如乘法的成本遠大於計算加減法的成本，此演演算法的速度會比原來的演演算法快。在現代先進的電腦上，乘法和加法花的時間可能相同，那麼此演演算法就沒有速度上的優勢。若使用浮點數計算，此演演算法可能會有精準度下降的問題，因此需要取捨。

在FFT（或是其他的線性映射），若是乘以固定係數（c + di，在FFT中稱為旋轉因子）的複數乘法，其中二個加法（d−c和c+d）可以事先計算，需要即時計算只剩三個乘法以及三個加法。不過若是配合有浮點運算器的處理器，加法的速度和乘法相當，因此減少乘法，增加加減法的演演算法，沒有速度上的優勢。

## 大數字的快速乘法演演算法

針對二個$n$位數數字的乘法，哪個乘法演演算法的速度最快？

### Karatsuba乘法

若需要計算到上千位數字相乘的系統，例如電腦代數系統或高精度計算函式庫，長乘法的速度太慢。因此會使用Karatsuba演算法，此乘法演演算法是Anatoly Karatsuba在1960年發現的，在1962年出版。此乘法演演算法的精神在於觀察到二位數的乘法可以不用傳統四個乘法的作法，可以只用三個乘法完成。這是分治法的例子。假設要將二個二位數m進位數字 x1 m + x2和y1 m + y2相乘：

- 計算x1 · y1，稱此結果是F
- 計算x2 · y2，稱此結果是G
- 計算(x1 + x2) · (y1 + y2)，稱此結果是H
- 計算H − F − G，稱此結果是K，此數字等於x1 · y2 + x2 · y1
- 計算F · m2 + K · m + G.

若要計算三個m進位數字的乘積，又可以用上述的技巧，使用遞歸的方式進行（但需要改為其他的進位），此方式。只要計算出乘積數字，再將數字相加，需要n個步驟。

Karatsuba演算法的時間複雜度是O(nlog23) ≈ O(n1.585)，此方式明顯的比長乘法要快。因為遞歸產生的額外計算，若n較小時，Karatsuba算法會比長乘法慢，因此在n較小時，會切換回長乘法進行計算。

Karatsuba演算法是第一個時間複雜度漸進的比長乘法快的演演算法，也可以視為是快速乘法理論的基礎。

1963年時，Peter Ungar建議將m改為i，以產生快速複數乘法的演演算法。若要計算 (a + b i) · (c + d i)，可依以下步驟進行：

- 計算b · d，稱結果為F
- 計算a · c，稱結果為G
- 計算(a + b) · (c + d)，稱結果為H
- 結果的實部是 K = H − F − G = a · d + b · c
- 結果的虛部是 G − F = a · c − b · d

如同上述複數乘法演演算法所述的一樣，需要三個乘法以及五個加減法。

### Toom-Cook 演算法

另一種乘法的演演算法是Toom-Cook演算法，或稱為Toom-3演算法。Toom-Cook演算法將每一個要相乘的數字切成幾部份，是Karatsuba演算法的推廣。三路的Toom-Cook演演算法可以針對大小為$3N$的被乘數和乘數進行乘法，其成本是大小為$N$的被乘數和乘數乘法的5倍，因為運算加速的係數是$9/5$，Karatsuba演算法的加速係數是$4/3$。

若用遞迴的方式將數字分成更多份，可以再縮減計算時間，但位數管理以及加法的成本也會增加。若位數到上千位數，一般會使用傅立葉轉換，速度多半會比較快，若位數更多，傅立葉轉換的優勢更加明顯。

### 傅立葉變換法

<image>

Strassen（1968年）曾提出用快速多項式乘法來計算快速整數乘法的基本概念。後來在1971年由Strassen和Arnold Schönhage找到理論和實務的確認，所得的就是Schönhage-Strassen演演算法。

## 下限

若用單一處理器將二個n進位的數字相乘，時間複雜度有一個trivial的下界Ω(n)，沒有更接近實際應用的下界。乘法在ACC0以外（p為任意整數），意思是不存在固定深度，多項式（甚至是次指數）大小，以AND、OR、NOT、MODp電路組合成的電路可以計算乘法。

## 多項式乘法

上述的乘法演演算法也可以用來計算多項式的乘法。例如Strassen演演算法就可以用來計算多項式的乘積。而Kronecker替代也可以將多項式的乘法轉換為二個整數的乘法。

以下是用長乘法來計算多項式乘法的例子：

```
 14ac - 3ab + 2 乘以 ac - ab + 1

 14ac  -3ab   2
   ac   -ab   1
 ————————————————————
 14a<sup>2</sup>c<sup>2</sup>  -3a<sup>2</sup>bc   2ac
        -14a<sup>2</sup>bc         3 a<sup>2</sup>b<sup>2</sup>  -2ab
                 14ac           -3ab   2
 ———————————————————————————————————————
 14a<sup>2</sup>c<sup>2</sup> -17a<sup>2</sup>bc   16ac  3a<sup>2</sup>b<sup>2</sup>    -5ab  +2
 `=======================================`<ref></ref>
```

## 外部連結

### 基本運算

- The Many Ways of Arithmetic in UCSMP Everyday Mathematics
- A Powerpoint presentation about ancient mathematics
- Lattice Multiplication Flash Video

### 進階演演算法

- Multiplication Algorithms used by GMP