高速乗算法の設計と実装(3)

梅谷 武
作成：2001-09-20 更新：2005-04-20
メルセンヌ素数を標数とする２次体上の離散フーリエ変換を使った32bit算術演算器向け高速乗算法を設計し、Pentiumへの実装を試みる。高速乗算法(2)においては係数を分解したが、高速乗算法(3)は係数を分解しないで計算する。
IMS:20010920001; NDC:412.1; keywords:離散Fourier変換,メルセンヌ素数,２次体;
目　　次
1. 高速乗算法(3)
1.1 高速乗算法(3)の構造
参考文献
1. 高速乗算法(3)
1.1 高速乗算法(3)の構造
Step 1. 数のP進表現
　基数をP=232とし、長さを
N=2n,0＜n≦25
とする。正の整数a,bが次のようにP進表現されるものとする。
a
＝
aN-1PN-1+…+a1P+a0,0≦ai＜P
b
＝
bN-1PN-1+…+b1P+b0,0≦bi＜P
この積について0≦ab＜PNが成り立つと仮定する。
Step 2. 多項式としての積の計算
　PK≡1(mod PN-1)より、
c＝ab≡N-1
Σ
r=0
(
Σ
s+t≡r(mod N)		asbt)Pr
(mod PN-1)
となる。この係数cr=∑s+t≡r(mod N)asbtの大きさを評価すると、
0≦cr≦N(232-1)2＜289-1
となる。
Step 3. 離散フーリエ変換による係数の計算
　crは、(N, ZM89[i], ζ), ζ=(6+i)t, t=290-n(288-1)型の離散Fourier変換を利用して次のように計算する。
F(a)k
=
N-1
Σ
s=0		asζsk　　(inZM89[i])
F(b)k
=
N-1
Σ
t=0		btζtk　　(inZM89[i])
cr
=
N-1N-1
Σ
k=0		F(a)kF(b)kζ-kr　　(inZM89[i])
Step 4. 桁上げ処理
　最後にcr, r=0,…,K-1に桁上げ処理を施すことで、c=abのP進表現が得られる。
参考文献
計算数論
[N1] 和田 秀男, コンピュータと素因子分解, 遊星社, 1999
[N2] 和田 秀男, 高速乗算法と素数判定法, 上智大学数学教室, 1983
[N3] 梅谷 武, 離散Fourier変換
[N4] 梅谷 武, ストラッセン-ショーンハーゲ法
[N5] Daniel J. Bernstein, Multidigit multiplication for mathematicians
[N6] 梅谷 武, 高速乗算法の設計と実装(1)
[N7] I. S. Reed, T. K. Truong, The use of finite fields to compute convolution, IEEE Trans.IT-21, 208-213, 1975
[N8] I. S. Reed, T. K. Truong, Complex integer convolutions over a direct sum of Galois fields, IEEE Trans.IT-21, 657-661, 1975
[N9] 梅谷 武, 高速乗算法の設計と実装(2)
算法
[S1] 野下 浩平, 高岡 忠雄, 町田 元, 基本的算法, 岩波書店, 1983
[S2] D. E. Knuth(中川 圭介訳), 準数値算法―算術演算, サイエンス社, 1986
[S3] H. J. Nussbaumer(佐川雅彦他訳), 高速フーリエ変換のアルゴリズム, 科学技術出版社, 1989
[S4] David H. Bailey, The computation of π to 29,360,000 decimal digits using Borweins' quartically convergent algorithm
[S5] Mikko Tommila, Number theoretic transforms
ディジタル信号処理
[D1] 電子情報通信学会, ディジタル信号処理の基礎, コロナ社, 1988
[D2] G. A. Jullien, Number theoretic techniques in digital signal processing
[D3] G. A. Jullien, Residue arithmetic with application in digital signal processing