暗号計算機屋のブログ

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

ICF3-Fのモンゴメリ乗算器が高速である仕組み

ICF3-Fは1999年のICF3FPGAに実装するために設計を、はじめたバージョンです。

MITの研究者が公開鍵暗号の高速化を可能にする乗算器を発明したとか、そういう話題を、たくさん輸入しはじめていると感じたので、僕が海外のパクリではないことを証明するために、僕のICF3-Fの種を明かすことにします。 ただ成功すれば商品化をしようかと考えているので、その辺りは考えてください。

基本原理は1999年のICF3と同じです。ICF3は、この教科書の「14.36 アルゴリズム モンゴメリ乗算」を実装しています。サンプルとしてFreeでダウンロード可能なので確認してみてください。

cacr.uwaterloo.ca

このアルゴリズムは基数(b)がパラメータになっています。1999年のICF3では基数をb=2としています。ICF3-Fでは、これをFPGAの持つ乗算器に合わせます。XilinxのDSP48E1は18x25bitの乗算器なので基数 b=218にします。基数が小さいうちは、小さい演算ブロックをいっぱい並べるような感じでICF3のb=2と同じように実装ができます。 ただしb=2よりも大きい基数を使うと逆数を演算する必要がでてくるため、少し複雑になります。ICF3-FではDSP48E1を1個利用して逆数専用の演算器を用意することになりそうです。またuiの生成も、それなりに複雑になるのでDSP48E1を1個利用してui生成専用の演算器を用意します。RSA 2048bitではCRTを使うことで1024bitのモンゴメリ演算器で演算が可能です。1024bitのモンゴメリ乗算器を24bit間隔で分割して小さいブロックにします。 1ブロックは1個のDSP48E1で構成されます。

ICF3では1024bitのモンゴメリ乗算器を2個装備するので使用するDSP48E1の個数は

( 1 + 1 + (1024÷24 ) ) × 2 = 90個

つまりXilinxFPGA秋葉原秋月電子通商で1個 9,800円 (2018年5月8日現在)で販売されているCmod A7-35Tに実装できるかもしれないということです。

14.36のアルゴリズムの話に戻します。恐らく14.36の2.2の式

A ← ( A + xi・y + ui ・m) / b

Aは約1024bit長なので、この加算キャリーに時間がかかって性能がでないと思った人が多かったのかもしれません。過去のIEEEの論文では、ASICの乗算器と同じような実装で、加算キャリーの問題を回避しているものが多いです。ポイントはAは加算されていきますが途中で参照されるのはAの最下位の桁だけということです。つまり基数b=218であれば、Aの最下位の18bit分だけ計算が完了していればいい。

余談になりますが、1999年のICF3が優秀なのは、この加算キャリーの対策の他、3の式

If A ≧ m then A ← A - m

を、非常に少ないゲート数で実装できたことなどの理由があります。そしてCRTや楕円暗号が1999年の製品で、動作するということです。

モンゴメリ乗算をFPGAに実装する場合、ASICと全く同じ方法では、性能はともかく、リソースの消費が大きくなると、予想されるので、式2.2のWallace Treeで実装する部分をブロックごとに全加算します。するとブロックは24bit幅ですが、18bit×24bitの乗算があるので、ブロックの演算結果は44bitになる。18bit右シフトして、ブロックから上位にはみ出す2bitを上位ブロックに送り、下位にはみ出す18bitを下位に送る。当ブロックでは、上位からはみ出してきた18bitと、下位からはみ出してきた2bitをマージして24bitのレジスタに入れる。次のサイクルで当ブロックの演算結果24bitと加算するとブロック内のAは25bit 1レジスタとなる。

もう一度、説明すると25bit + 18×24bit + 18×24bit = 44bitになり、延々、ブロック内のAは25bitで足りる。

この結果、加算キャリーは上位1ブロックまでで1024bitの加算キャリーを毎サイクル必要とすることはなくなり、並列に高速に演算できるようになるのです。 (2018年5月14日 公開初日の説明が、少し甘かったので修正しましたが、結論に違いはないように思います。 ASICだと自由に設計できるので解がユニークに決まってくるのですが、FPGAだとDSPに合わせる都合、いろいろになって少し戸惑いました。ブロック毎に全加算して楽をしようと思ったところとか。DSP48E1だと17bitシフタがあるからそれを有効に利用するなら218よりも217になるし、変数Aを48bitのCレジスタに割り当てることができれば25bit幅でも良くなる。)

結論

ICF3-Fは、ほぼ1999年に発明されたいた。真に難しいのは演算器を並列動作させることよりも、演算の全工程が可能である実装、暗号プロセッサの発明が重要であること。 FPGAに実装してみたら、かなり効率的に実装できそうなので、ASIC並みの性能が期待できる。FPGAは需要に応じてRSA 4096bit対応させてたり、楕円対応で256bitに変更できるので、FPGAのハードウェアを効率的に運用できるのでASICを超える可能性すらある。ASICではRSA 4096bitと楕円256bitを同時にできないといけないので無駄が生じる。

インターネットは、昨年から常時httpsが必須となりSSL演算処理のコストが増加したが、Cmod A7-35TをUSBに刺すだけで、処理能力が向上するため、性能向上分のWebサーバーのコストが浮く。常時httpsで使われるopensslが非同期処理に対応しているので現在販売されているUSB型のFPGAでいけそう。つまり

浮いたサーバーコスト >  Cmod7-35T を含めた販売価格

であれば商売が可能だということです。FPGAが実際の生活に役に立つあたりも、いいと思っている点です。

FPGAを使ったSSLアクセラレータのICF3-F

これから設計して性能が出そうなら売れるものになるのかも、という話で具体的な計画があるわけではないです。 1999年の暗号LSI ICF3をベースとしてFPGAに特化したSSLアクセラレータを考えていきます。これまでFPGASSLアクセラレータでは性能が出せず商売にならないと思っていたのですが、FPGADSPを多数もっていることを知り考えが変わりました。

演算器の性能でいえばGPUのほうがはるかに上なのですが複数のSSLの演算を同期しなければならず利用ケースは限られます。FPGAではモンゴメリ乗算をうまく使うことで多数のDSPを効率的に使用できるはずなので小規模なサーバーでも効果があることを見込める。

これまで、そういう研究は10年以上前からありましたが商用化に成功しているものは、あまり見当たらない。あまり探していないので、もうあるかもしれないですが。商用化を困難にしている原因のひとつは、SSLの演算の全てをFPGAだけで処理することが難しいからではないかと「勝手に」思います。

そこでICF3をベースにFPGAに特化したSSLアクセラレータを作れば、FPGAだけで処理できるようになる。FPGAだけで処理できるようになるとUBS接続の安価なFPGAでも小規模サーバーのSSLに貢献できるようになり商業化が可能になるかもということです。

実際にやってみて性能が出れば。ですが。(^^;; うーん、ICF3-Fとして検討してみるか。

2018年4月23日 19:35追記

実際にやってみないと、わからないところはある。技術的な説明をすると、ICF3は基数2で幻のICF4は高基数に振り切っている。FPGAではDSPのビット長にあわせて基数を2から、もう少し大きくする。基数2にない面倒なことが起きるかもしれないが、そのくらいならICF3のように全体を閉じることができないかと考えている。 これが可能になっているのは僕が、ハード、ソフト、暗号を1人でできるからなのだと思う。仕様変更も頭の中で3秒くらいで完了するし。

IT業界で数学の才能が何の役に立つか

取り合えず僕に数学の特技があると思って聞いてもらえるといいと思います。😉

1999年の暗号LSI ICF3の開発で僕がモンゴメリ乗算を採用することに反対した理由を説明すれば、数学の才能が役に立つということが、わかるように思います。

ICF3を開発した部署は、ハードウェアの開発部なので電気・電子系の学科の人が多いのですが、現在の教育のようにRSAの数学について、知っている人は、いませんでした。研究所とは異なり頭を空っぽにして、全力で、最小限の工数で、開発を達成させることが仕事なのです。特に隣の大型コンピュータの開発部隊と異なり1年に1機種を開発するという無謀なことをやっていました。隣では4年に1機種というペースでやっています。なので勉強なんてさせてもらえず成果だけ出すという苦しいことをやらされます。1年に1機種が可能なのは、それなりの学歴の年配の方々が控えていて、できなければ、切り捨てられてしまうからです。

そんな中、新規にRSA暗号LSIを開発することになった。検討を進めて方式が決まりそうになったところで日立のシステム開発研究所から海外の暗号の教科書がFAXで送信されてきた。22ページだけ抜き取ったものでした。モンゴメリ乗算を使ってRSAを実装できないかということだった。研究所の人が必要だと思った、バラバラの箇所の22ページで、最後の2ページにモンゴメリ乗算の式が、書かれていた。「14.3.2 Montgomery reduction」のところです。 そこから2日かけてRSA暗号の計算方法を、作り上げたわけです。実は教科書のアルゴリズム14.94にRSA暗号の計算方法があったのですがFAXには入っていませんでした。研究所の人は、ちゃんと読んで送っているのか、疑いました。😀 ただし、ハードウェアでは14.94の方法では難しくて、どのみち自身でRSA暗号の計算方法を作り上げる必要がありました。

モンゴメリ乗算の公式14.36 --- x・y・R-1 mod m

R-1は、引く1ではなくてマイナス1乗です。RSA暗号の計算で、R-1が消えるように、計算していけば、答えがでるだろうと予測して、C言語で、試してみたところ、正解が得られました。

ここで数学の才能の有無によって、判断に違いが出てくるのです。R-1という記号が、自分が知っている-1なのか、あらゆる値で、正しい答えが計算できるのか、教科書の読んでないところに、特異点が存在するみたいなことが書かれてないか。

大学1年のときの数学で整数論を習っていたのを思い出した。10年前に習ったことで、加減、乗除を定義して、それが閉じるのか、どうか、という講義だった。モンゴメリ乗算という定義が、閉じているのか、それを確認する必要があった。

ただコンピュータを作りたいと思う電気工学科が、整数論を、しっかり覚えていたわけではなく、特異点が存在しないか、おろおろ、やっていた。最初、感覚的に問題を感じ、自分がわかっていないことが、わかったのだが、徐々にわかっていないところが、何なのか、明らかになり整数論にたどり着いた。

ソフトウエアの開発では、間違いがあれば、バージョンアップをすればいいという考え方もありますが、ハードウェアの開発では、1回の間違いが致命的な損害になることもあり、数学の才能はIT業界で役に立つということになると思います。

以下はFAXで送信された教科書です。サンプル無料とありますが、全部あるようです。

cacr.uwaterloo.ca

2018年4月19日追記

この話は「車輪の再発明」ではなくて、なんかやってみたら正しい答えになるから大丈夫だと勘違いをしている可能性があることを問題として認識できるか?という点です。 自身で計算方法を作っていますから、公式が、正しい条件で運用されているかを、厳格に考えなければならない。 割算(÷)の記号が使われていても、実際の演算はユークリッド互除法だったり、べき乗剰余演算だったりするのです。大学の数学は。

ついでの話になりますが、この国には自分で考えたほうが、海外よりも、できる回答になる人がいることを、頭の中に入れてもらえないかと思うのです。

1997年 東大の「疑似非同期式マイクロプロセッサ」

僕は1994年に日立 中央研究所 超高速プロセッサ部に入ったが、次の年には超高速プロセッサ部は無くなった。プロセッサを作るやつは「犯罪者」であり受刑者のような雰囲気だったのは確かだ。そして僕はIBMのCPUを買うための仕事にとりかかった。

1996年、大学などの研究機関向けにVDECが設立されLSIの試作、研究が行えるようになったようです。VDECの年報、1996年と1997年を読みました。東大のプロセッサのための研究が、いくつか掲載されています。 この頃は、パイプライン段数を増やす、スーパーパイプラインという思想で、周波数を上げる技術の研究が盛んだったようです。加算器をパイプライン化する研究や、加算する値によって計算時間が異なることを利用し、1~3サイクルで加算をする「疑似非同期」というアイディアがあったようです。

僕のICF3-Vの「疑似パイプライン」は面積当たりの性能を向上させるためパイプライン段数を削減しつつ、パイプライン効果を得る方法なのだが、それは、おいておいて、実は「疑似非同期」も検討していた。

東大の疑似非同期の対象は加算器で、入力値によって変化する終了検出が必要だった。僕の疑似非同期は、除算が対象。除算支援ハードのために全体の周波数が落ちているのだ。除算演算をする場合、32bitの2命令を1サイクルとするみたいな、決まり事を作るだけ。値による終了検出はしないので、簡単。除算の性能が下がるが、決まり事を作るだけでハードウェアの変更の必要はない。同一のハードで、除算の性能が下がっても、周波数を上げたい場合、ソフトウェアのほうで対応する。

VDECの年報は基本的には教育を含め、学生が短期間に試作しているものだが、中には、興味あるものもあるかもしれませんね。 1998年以降も、ちょっとみてみようかな。

http://www.vdec.u-tokyo.ac.jp/REPORT/vdecreport.html

MicroBlaze性能評価

目的

手っ取り早く自作CPUの性能を評価するためMicroBlazeと比較できるようにする。ベンチマークは自作CPUの性能を測るプログラムを簡単にするため、あまり意味はなく四則演算や論理演算、シフト命令が適当に入るような小さいものを自作した。FPGAの経験は浅いので結果についての取り扱いは注意されたし。

MicroBlazeの構成

FPGAバイス : Artix-7 35T Arty

性能型 面積型 最小型 MCS型
Predefined
Configuration
microcontroller Minimum Area Minimum Area -
implementation
optimization
PERFORMANCE AREA AREA AREA
クロック 166MHz 166MHz 166MHz 166MHz
WNS 0.04 ns 0.01 ns 0.02 ns 0.01 ns
Debug Module
Enable
Barrel Shiffter
×
FPU × × × ×
Integer Multi MUL32 MUL32 × ×
Integer Divider × ×
Optimaization PERFORMANCE AREA AREA -
Enable Branch
Target Cache
× × × ×
S-RAM 32KB 32KB 32KB 32KB
追加IP AXI GPIO
AXI Uartlite
AXI GPIO
AXI Uartlite
AXI GPIO
AXI Uartlite
-

GPIOはArtix-7 35Tの基板上のLED 4個と接続されている

実装面積

性能型 面積型 最小型 MCS型
LUT 1742 1512 1401 985
FF 1582 1274 1142 747
BRAM 8 8 8 8
DSP 3 3 0 0
LUTRAM 140 218 216 203

測定結果

パラメータ 性能型 面積型 最小型 MCS型
ループ数10000000
除算割合 10
7.9 秒 10.5秒 41.7 秒 51.0 秒

ベンチマーク プログラム

#include <stdint.h>
#include "xparameters.h"
#include "xil_printf.h"

#pragma GCC optimize("O2")

uint32_t xor(void) {
  static uint32_t y = 2463534242;
  y = y ^ (y << 13); y = y ^ (y >> 17);
  return y = y ^ (y << 15);
}

uint32_t xors(int loop,int div) {
    int i;
    uint32_t a,b,c;
    uint64_t w;
    uint32_t r;

    r = 0;
    for(i=0 ; i<loop ; i++) {
        do {
            a = xor();
        } while(a==0);
        do {
            b = xor();
        } while(b==0);
        do {
            c = xor();
        } while(c==0);

        w = a*b;
        if(i%div==0) w /= c;
        r += w;
    }
    return r;
}

#pragma GCC optimize("O0")

int main() {

    xil_printf("START!\n");
    uint32_t r = xors(10000000,10);
    xil_printf("%08X\n",r);

    return 0;
}

ベンチマークの演算結果

gccのオプションでハードウェア乗算器、除算器、バレルシフタの有無の設定ができる。MicroBlazeの構成をもとに正しくgccのオプションを設定する。乗算器がない構成で乗算命令を行うと演算結果が0になる。正しい演算結果になっているか確認する。測定したパラメータでは 次のような値になる。 0x3EEF270E

乗算の部分は32bit×32bit = 64bit のつもりで作成したが64bitにキャストしていなかったため32bit × 32bit = 32bitというようにCコンパイラは解釈するみたいです。結果の性能は32bitのままデータの再収取はしませんが、64bitにキャストした場合の期待値は

0xFC3BB111

Intel CPUとの比較

Intel Celeron E7500(2.93GHz)では0.13秒

2930 ÷ 166.667 ≒ 17.58倍 E7500の周波数を166.667MHzに換算した場合は2.3秒

考察

乗算器、除算器がある場合、性能は格段に向上する。乗算器のゲート数はASICでは大きなものになるはずだがFPGAではDSPで実装されるため注意する。性能型はパイプライン段数5で面積型はパイプライン段数3なので、性能型のほうが周波数が高くてもいいはずなのだが、今回は同じ166MHzだった。Xilinxのサイト見ても周波数の差は少ないので合っているものとも思われる。乗算命令とシフト命令のレイテンシが性能型と面積型で異なることが、今回の測定結果に反映されたのかもしれない。乗算命令のレイテンシは性能型で1、面積型で3。シフト命令のレイテンシは性能型で1、面積型で2。 MCS型はMicroBlaze MCSのこと。機能固定で面積をさらに小さくしたものと思われる。機能面では最小型とバレルシフタの有無が異なる。この結果、同じ周波数だが、バレルシフタがある最小型のほうが高性能になっているものと思われる。

追記 2018年4月14日 次のURLでMicroBlazeの乗算性能についての説明があって参考になります。

qiita.com