Eigen ー C++で線形代数を!

fukurou
  • ????????????????????

このエントリーおよび続くエントリーで、C++ で書かれた線形代数ライブラリ Eigen を使うためのメモなどを書いてみようと思います(このエントリーは旧ブログからの転載です)。

とはいっても本家のほうに非常にわかりやすいマニュアルが整備されているので、「この部分はマニュアルを見てね」的な部分も出てくるかもしれませんが。
最終的には微分方程式とか最適化とかの応用問題について書いてみたいと思います。

(現時点(2011/2)で私が触っているのはバージョン 2.0 ですが、もうすぐバージョン 3.0 が出るようです)

なぜ Eigen?

私は数値計算のプロではなく、ユーザー側の人間なので以下のような点に注目していろいろと線形代数ライブラリを漁ってみました。

  • 自分の頭の中や論文の数式をそのままコーディングできるような直感的でシンプルな API
  • 世界最高でなくてもよいが、比較的優れたパフォーマンス
  • Fortran でなく(できれば C でもなく)C++ での実装
  • 疎行列 (sparse matrix) のサポート
  • 割と応用的な操作(連立方程式、固有値問題、特異値分解など)のサポート
  • ある程度の世の中での使用実績
  • BSD or LGPL ライクなライセンス

C++ での線形代数ライブラリの実装には

  • uBLAS (boost)
  • Eigen
  • Blitz++
  • CPPLAPACK
  • LAPACK++

などがあるようです。

こうしてライブラリを眺めてみると前者三つは「式テンプレート型」、後者二つは「BLAS/LAPACK ラッパー型」ですね。

これらのライブラリと上のような個人的評価基準をもとに、頭の中でゴニョゴニョしたり、テストプログラムを走らせたりした結果、Eigen を使うことにしました。

Eigen の特徴

以下、本家から引用してみます。上であげた個人的に注目するポイントはほぼクリアしていることがわかります。

シンプルな API

The API is extremely clean and expressive, thanks to expression templates. Implementing an algorithm on top of Eigen feels like just copying pseudocode. You can use complex expressions and still rely on Eigen to produce optimized code: there is no need for you to manually decompose expressions into small steps.

(意訳)式テンプレートのおかげで非常にシンプルな API を実現できました。Eigen を使えばほとんど擬似コードをコピペする感じでアルゴリズムを実装できます。複雑な式も Eigen が自動的に最適化してくれるため、ユーザーは式を(たとえばテンポラリーな変数を減らすために行う)小さなステップに分割する必要がありません。

(コメント)たとえば z=Ax+y みたいな式は、式テンプレートを使わないと行列とベクトルの積を計算して、それをテンポラリーな変数に代入してそれと y の和をとる、みたいな無駄が発生します。式テンプレートの使用は右辺を(行列×ベクトル+ベクトル)型、として扱うことでテンポラリーな変数を削減し、実行時間を短縮する、というイメージでしょうか?

速度

Expression templates allow to intelligently remove temporaries and enable lazy evaluation, when that is appropriate — Eigen takes care of this automatically and handles aliasing too in most cases.
Explicit vectorization is performed for the SSE (2 and later), AltiVec, and ARM NEON (new in Eigen 3) instruction sets, with graceful fallback to non-vectorized code. Expression templates allow to perform these optimizations globally for whole expressions.
With fixed-size objects, dynamic memory allocation is avoided, and the loops are unrolled when that makes sense.
For large matrices, special attention is paid to cache-friendliness.

(意訳) 式テンプレートはスマートに一時変数を取り除き、遅延評価を行うことを可能にします。SSE, AltiVec, そして ARM NEON (Eigen 3 以降) などの各種ベクトル演算命令セットに対して明示的にベクトル化がなされます。式テンプレートはこれらの最適化を式全体に対してグローバルに実行することを可能にします。固定サイズの行列・ベクトルに対してはループが展開されます。巨大な行列に対してはキャッシュフレンドリーになるように特別な配慮がされています。

(コメント) 式テンプレートはすばらしい技術だと思いますが、あのコンパイルエラー時のメッセージはなんとかならないものでしょうか・・・

スパース行列のサポート

both dense and sparse (the latter is still experimental) matrices and vectors.

(意訳) 密行列と疎行列両方のサポート(疎行列に関しては実験的)

(コメント) 現時点で実装されている疎行列の格納方法は compressed matrix のようです。疎行列にはさまざまな要素の格納方法があります。たとえば boost::ublas では連想配列を使った疎行列 (mapped matrix) も実装されていたりします。

実績

科学技術計算、CG、ロボット、ゲームなどいろいろと使われているようです。有名どころでは Google が機械学習に使っている、という実績もあるそうです。

ライセンス

GPL/LGPL の選択式。

サポートしている操作

連立一次方程式(LU分解のみでクリロフ部分空間法的な反復法はサポートなし)、特異値分解、QR分解、コレスキー分解。

次回予告

次回以降は具体的に Eigen を使ったサンプルなどを書いていく予定です。

はてなブックマーク - Eigen ー C++で線形代数を!
Pocket