はじめに
様々な用途で大規模行列の固有値計算をする必要に迫られており、簡単に使えるライブラリを調べてみた結果、やっぱり課金しただけあってIntel MKLが一番だなという(ある意味)当然の結論に落ち着きました。何度も同じことで悩みそうなので備忘録を書いておきます。殴り書きに近いのでそのうち修正します。
解きたい問題
行列の低ランク近似がやりたいだけです。
任意の実数対称行列$\mathbf{A}\in\mathbb{R}^{N\times N}$について、
$\mathbf{A}\mathbf{t}=\lambda\mathbf{t}$
を満たす固有値・固有ベクトルのペアを、固有値の絶対値が大きい順に$m\le N$個だけ取得したい。
使えそうな方法
べき乗法(power method)・subspace法
たぶん古典的にはsubspace法が一番有力で、70年代に作られた割に振動解析の分野など疎行列の固有値計算の世界ではまだまだ現役を誇る。でもsubspace法のライブラリが(調べた限り)不思議と見つからなかったので検討外。
Arnodi法
超有名な手法でググると解説がたくさんでてくる。丸め誤差に弱いなどの問題があり、逐次計算の中で何度かリスタートなどの処理をする必要がある。implicit restart Arnordi法は最大の成功例であり、このアルゴリズムを使ったFORTRANライブラリARPACKがよく使われている。ただ、90年代にできたプログラムなので今の環境では工夫しないとコンパイルが通らずインストールできない。C++のラッパーも公式が提供しているが、こいつもいかんせん古くて動かない。平成も終わるこの時期にわざわざこれに頼らなくてもいいかなと思う。
ARPACK公式
https://www.caam.rice.edu/software/ARPACK/
実対称行列の固有値問題に対するARPACKの適用事例(日本語)
http://blog.beam2d.net/2011/12/eigenarpack.html
上記著者の方によるC++ラッパー
https://github.com/beam2d/arpaca
FEAST
正直調べるまで存在自体を知らなかった。全く内容を理解していないが、後述するMKLの解説の中では量子力学の知見から着想を得て生まれたアルゴリズムとある。何を言っているのかわからない。実際、最初にこのアルゴリズムが提案されたのはPhysical Review Bなので、おそらく本当。FEAST LIBRARYとしてオープンソースのライブラリが2009年に登場し、この移植版がIntel MKLに入っている。exampleの場所がわかりにくいが、
MKL/examples/solvers_eec/
の中に、疎行列・密行列・複素行列のそれぞれの場合での固有値計算の例がある。
FEAST公式
http://www.ecs.umass.edu/~polizzi/feast/
Intel MKL公式
https://software.intel.com/en-us/articles/introduction-to-the-intel-mkl-extended-eigensolver