本当は怖いHPC

HPC屋の趣味&実益ブログ

固有値ソルバの初歩の初歩、ヤコビ法を書いてみた

行列の固有値を求めるアルゴリズムとして一番初歩的なのが、ヤコビ法です。実対称行列の固有値を求めるためのアルゴリズムで、実対称行列の固有値は実数であることから実装も簡単です。

理屈としては、

  1. 対角行列の固有値は対角成分そのもの
  2. 直交行列 Pを用いて \Lambda = P^{-1}AP と書けるとき、A\Lambda固有値は一致する

という2つの定理を使って、適当な直交行列を準備して次々とAを変形していって、最終的に対角行列にできれば固有値が求まる、という寸法です。詳細は数値解析入門 を参考にしました。三角関数を使って直交行列を巧妙に作ることによってAの非対角要素を次々と0に消し去っていくのが鮮やかですね。

なお、Pが直交行列であるとは、 Pの転置行列をP^{T}と書くとき、PP^{T} = E が成り立つ行列のことです。

ソースコードはこちらhttps://github.com/keisukefukuda/numerical/tree/master/eigen_jacobi

実験してみた結果、行列のサイズが大きくなると収束までの繰り返し回数が非常に大きくなってしまうので、巨大な行列に対してはもっと高度なアルゴリズムを使う必要があります。

次はQR法もしくはHouseholder-Givens法を勉強してみます。

【広告】