“Fast matrix multiplication is stable,” but…

 今年(2014年)の3月に応用数理学会連合研究会で多倍長Strassenアルゴリズムについて喋った時に,”Group-theoreticアルゴリズム”の話がちょろっと出たので,Strassenの件が落ち着いた時にちょっと調べてみたのである。・・・なるほど確かに「一時は」盛り上がっていたようだが,実装に絡む大事な論文が消えてしまって以来,どーも最近雲行きが怪しいようなのである。とはいえ,面白そうな話題ではあるし,ちょっとずつ勉強していくには良さそうなものなので,コメントいただいた方に対して「一応は調べてみましたよ」という報告がてら,この記事を起こしたという次第である。「実装に絡む大事な論文」のその後や”Group-Theoreticアルゴリズム”の動向をご存知の方がいらっしゃったら情報も頂けると幸いである。

 つーことで多倍長Strassenの論文掲載が決まったので,ソースも外に出した。ろくすっぽ整理してない代物なのだが,ちょろっと引き合いもあったこともあって,ちゃんと動いているという証拠物件なので早めに出すことにしたのである。
 ちなみにこれ,比較対象のため倍精度の行列積もIntel Math Kernel(IMKL)提供dgemm関数を利用して実装しており,その結果については紀要原稿に書いてある。ま,メインは多倍長行列積の方なので,あくまで比較用なのだが,存外実装が難しいなぁという感想を持った。特に2^nでない行数・列数の時の処理が大変で,今のところは多倍長に合わせてパディング(0要素を詰め込んで2^nの倍数の行数・列数にする)とピーリング(2^nの倍数をはみ出る部分はブロック化して別に処理する)処理を組み合わせて実装しているが,ATLASベースのStrassen実装報告にある通り,倍精度ぐらいだとパディングだけでサイズ合わせをした方が無難かもしれない。ちなみに奥村先生の「C言語によるアルゴリズム事典」「Javaによるアルゴリズム事典」にもStrassenのアルゴリズムの話題があるが,1988年のBailey論文での成果を元にしている「2048×2048の行列積で5割高速」という記述をそのまま真に受けるとバカを見る,ということはATLASベースの実装とワシのやったIMKLベースの実装,両方の結果を見てもらえばご理解いただけるだろう。現在までのStrassenアルゴリズムのサーベイについてはワシの紀要原稿にもちらと書いたが,スポック博士の分厚い本の23章が一番コンパクトにまとまっているので参照されたい。

 つまり倍精度計算だと,キャッシュヒット率・SIMD命令の活用・メモリ操作関数のパフォーマンスによる影響が大きく処理時間に作用することになるので,演算量を多少減らしたところで,Strassenアルゴリズムを実行するために必要となる前処理・後処理をも含めると,全体として処理時間が減らせるかどうかはやってみないとよく分からないところが大きいのである。その点,元々四則演算の処理速度がトロい多倍長計算だと,演算量の減少がダイレクトに処理時間の減少に繋がるので大変分かりやすい結果になる。つーことで任意サイズの行列積の計算がStrassenアルゴリズムで可能になれば,LU分解の処理時間も最大3割減となるのである。逆に言えば,今時倍精度計算でLU分解に適用しても,PC単体のオンメモリで収まるサイズの問題程度では殆ど効果はないのでは,と思えるのである。Bailey先生のCray-2と今のx86環境とでは,キャッシュやSIMD命令をサポートしたアーキテクチャも大分異なっているので,同一視してはいかんのである。

 そんだけ,今の計算環境では高速だが前処理・後処理の複雑なアルゴリズムの実装は職人芸が必要であり,ましてやStrassenやWinogradアルゴリズムよりも複雑怪奇な”Group-Theoretic アルゴリズム”って奴の実装は相当な困難が予想される。これについてはオリジナルがこれで,それに対して誤差解析の観点から”Stable”という判定を下した(イマイチ疑問があるのだが)論文がこれであるらしい。後者については査読論文になっているようなので(これ),内容的にどうこういう気はないのだが,問題は後者の奴にある参考文献[8]。どうやら実装した結果についての報告らしいのだが,検索しても見つからないし,同様に実装した結果について記述したものが一つも見つからないのである。どなたかチャレンジして失敗したって記述でもあるなら分かるんだが,どーも今のところちゃんと実装した人が一人もいないんじゃないのかとワシは大いに疑っているのである。参考文献[8]が見つからないのは実装に失敗したか,作ってみたけど倍精度程度では高速にできなかったんじゃないのかな~,と今のところ考えている。ま,Strassenアルゴリズムの結果を見る限り無理もないかなぁとは思うのである。前処理と後処理にめんどくさそうなFFT?を使っていようだから,かなり難しいんじゃないのかしらん? まぁだからこそ「多倍長でやってみたらいいのでは」というご意見も分かるんだけど,さーて・・・?

 と,この先はワシが考えていけば済むことなので,少し時間をかけてアルゴリズムを解読してみるが,何分,かなり浮世離れした感のある研究テーマだし,その後の進展が聞かれないところを見ると,いろいろ実装上の困難もありそうである。だからこそのチャレンジングなテーマとも言えるけど,そういうのはもっと頭のいい人に取り組んでもらえないかなぁと思っているのである。もしどなたか既に着手されてそれなりの成果が出ているようでしたら,お知らせ頂ければ幸いなのである。あ,黙っていてくれと言われれば黙ってますのでこっそりでも結構でございます,はい。