暗号計算機屋のブログ

なにか思いついたことを不定期に更新。

モンゴメリ乗算の累積加算における分割加算の証明

はじめに

今年の4月にFPGAをはじめたときにFPGAにDSP(乗算器)が多数あることを知った。 これでモンゴメリ乗算が高速化できるわけだがアルゴリズムをそのまま実装するとビット長の大きな累積加算が、ビット長とともに周波数が下がり性能がでない。 そこで分割して加算しても、正しい加算結果になるような方法を考えた。これを使ったSSLアクセラレータ ICF3-F(商用版)の実装を急いでいるところだが、証明もしておかないと、いけないかと思って、急いで公開することにしました。産業スパイによる海賊版を抑止する目的です。ちなみに他にも応用ができる方法ではないかと思います。 この証明にコメントしたい方が、いらしゃるのかわかりませんが、もしあれば公開日(2018年9月26日)の２日後、28日の午後12時からでお願いします。とりあえず、その記録は残すようにします。

問題

非負の整数の変数Aを2進数表現すると、最下位ビットからrビットづつ区切ってs個のブロックに分割して表現できる。

\displaystyle A = \sum_{k = 0}^{s-1} a_k \cdot 2^{r k}

a_k がrビットの場合はAはr×sビットの数を表現可能だが、 a_kが(r+1)ビットあってもr×sビットの数を表現できる。 最上位ビットを e_k、それ以外を f_kとすると a_k = e_k \cdot 2^r + f_kになるが

　\displaystyle A = \sum_{k = 0}^{s-1} e_k \cdot 2^{r (k+1)} + \sum_{k = 0}^{s-1} f_k \cdot 2^{r k}

として全加算してやればAは冗長性がなくなりr×s+2個のビットの数になる。 このような冗長性をもったAを考える。

A = (A + t \cdot B + u \cdot C ) \div 2^{n} --- (1)

r \geqq n + 2 --- (2)

t,uはnビットの数。B、Cはr×s-2ビットの数。 \displaystyle B = \sum_{k = 0}^{s-1} b_k \cdot 2^{r k}、 \displaystyle C = \sum_{k = 0}^{s-1} c_k \cdot 2^{r k}と表す。

式(1)の右辺を計算した結果、任意のA、B、C、t、uにおいて、 同じ冗長性、すなわち a_kが(r+1)ビットある左辺Aに正しく格納できることを証明せよ。

証明

式(1)右辺のカッコ内は

A + t \cdot B + u \cdot C

\displaystyle = \sum_{k = 0}^{s-1} a_k \cdot 2^{r k} + t \times \sum_{k = 0}^{s-1} b_k \cdot 2^{r k} + u \times \sum_{k = 0}^{s-1} c_k \cdot 2^{r k}

\displaystyle = \sum_{k = 0}^{s-1} a_k \cdot 2^{r k} + \sum_{k = 0}^{s-1} t \cdot b_k \cdot 2^{r k} + \sum_{k = 0}^{s-1} u \cdot c_k \cdot 2^{r k}

\displaystyle = \sum_{k = 0}^{s-1} ( a_k \cdot 2^{r k} + t \cdot b_k \cdot 2^{r k} + u \cdot c_k \cdot 2^{r k} )

\displaystyle = \sum_{k = 0}^{s-1} ( a_k + t \cdot b_k + u \cdot c_k ) \cdot 2^{r k}

ここで p_k = a_k + t \cdot b_k + u \cdot c_kとする。t,uがnビットであるから p_kは最大でもr+n+2ビット。 p_kの最上位ビットから2bitを px_k、最上位ビットから2ビット、最下位ビットからnビットを除いたrビットを py_k、 最下位ビットからnビットを pz_kとする。

p_k = px_k \cdot 2^{r+n} + py_k \cdot 2^{n} + pz_k (0 \leqq k \leqq s-1) --- (3)

(1)の右辺は

\displaystyle \sum_{k = 0}^{s-1} p_k \cdot 2^{r k} \div 2^{n}

(3)より

\displaystyle = \sum_{k = 0}^{s-1} ( px_k \cdot 2^{r+n} + py_k \cdot 2^{n} + pz_k ) \cdot 2^{r k} \div 2^{n}

k=0の項をシグマの外に出して整数にすると

\displaystyle = px_0 \cdot 2^{r} + py_0 + \sum_{k = 1}^{s-1} ( px_k \cdot 2^{r+n} + py_k \cdot 2^{n} + pz_k ) \cdot 2^{r k} \div 2^{n}

\displaystyle = px_0 \cdot 2^{r} + py_0 + \sum_{k = 1}^{s-1} ( px_k \cdot 2^{r} + py_k ) \cdot 2^{r k} + \sum_{k = 1}^{s-1} pz_k \cdot 2^{r k - n}

\displaystyle = px_0 \cdot 2^{r} + py_0 + \sum_{k = 1}^{s-1} ( px_k \cdot 2^{r} + py_k ) \cdot 2^{r k} + \sum_{k = 1}^{s-1} pz_k \cdot 2^{r(k-1) + r-n}

\displaystyle = py_0 + pz_1 \cdot 2^{r-n} + \sum_{k = 1}^{s-2} ( px_{k-1} + py_k + pz_{k+1} \cdot 2^{r-n} ) \cdot 2^{r k} + ( px_{s-2} + py_{s-1} ) \cdot 2^{r (s-1)} + px_{s-1} \cdot 2^{r s}

ただしB、Cがr+n-2ビットなので必ず px_{s-1} =0。

\displaystyle = py_0 + pz_1 \cdot 2^{r-n} + \sum_{k = 1}^{s-2} ( px_{k-1} + py_k + pz_{k+1} \cdot 2^{r-n} ) \cdot 2^{r k} + ( px_{s-2} + py_{s-1} ) \cdot 2^{r (s-1)}

ここで式(2)の条件からr-nは2以上。 w_k = px_{k-1} + pz_{k+1} \cdot 2^{r-n} (1 \leqq k \leqq s-2)とすると w_kはrビットの記憶素子を持つ1つの変数としてまとめられる。

\displaystyle = py_0 + pz_1 \cdot 2^{r-n} + \sum_{k = 1}^{s-2} ( py_k + w_k ) \cdot 2^{r k} + ( px_{s-2} + py_{s-1} ) \cdot 2^{r (s-1)}

これは(1)左辺のAとして次にようになる。

\displaystyle A = \sum_{k = 0}^{s-1} a_k \cdot 2^{r k}

a_k = \left\{ \begin{array}{} py_0 + pz_1 \cdot 2^{r-n} ( k = 0) \\ py_k + w_k ( 1 \leqq k \leqq s-2 ) \\ ( px_{s-2} + py_{s-1} ) ( k = s - 1) \end{array} \right.

py_k、 w_kともにrビットの数なので加算してもr+1ビットの数である。k=0、k=s-1も同様にr+1ビット。 式(1)の右辺を計算した結果、同じ冗長性、すなわち a_kが(r+1)ビットある左辺Aに正しく格納できる。

補足

実際の実装では最後に、全加算を行って冗長性をなくしてr×sビットのレジスタに格納します。 このとき桁あふれが発生する可能性がありますが、一般のCPUのレジスタと同じで桁あふれはアルゴリズム側の責任。 モンゴメリ乗算では最終的な計算結果はCの2倍以下になるので全加算後、r×sビットのレジスタに必ず格納できます。 計算過程においても、上記証明によって、桁あふれすることなく計算できていることがわかります。 この問題では積和の数が２個の場合ですがn個にした証明など一般化できるように思います。いづれまた。