Eigenで線形代数

Eigenは線形代数,すなわちベクトルや行列を扱うことに特化したC++言語向けのテンプレートライブラリです. このライブラリを利用すればC++上でベクトルや行列を処理することが出来るようになります. Eigenは現在メジャーバージョンが3なので,以後それと分かるようにEigen3と呼ぶことにします. Eigen3は以下に挙げる9つのモジュールから構成されています.

Core: 行列と配列のクラス,基本的な演算.
Geometry: 幾何学的変換.
LU: 逆行列,行列式,LU分解とソルバ.
Cholesky: LL^TおよびLDL^TへのCholesky分解とソルバ.
Householder: Householder変換.
SVD: 特異値分解(SVD)と最小二乗ソルバ.
QR: QR分解とソルバ.
Eigenvalues: 固有値,固有ベクトル分解.
Sparse: 疎行列およびその基本的な演算.

また,"Dense"というモジュールにはCore,Geometry,LU,Cholesky,SVD,QR,Eigenvaluesモジュールすべてが含まれており,"Eigen"というモジュールにはDenseに加えSparseモジュールも含まれています. 基本的な行列計算をしたいが,どのモジュールをインクルードすれば良いのか迷った場合には,Denseモジュールをインクルードすれば間違い無いワケです. Eigenのモジュールすべてを利用したくなったらEigenをインクルードすれば良いワケです. Householderモジュールだけどこかに行ってしまった気はしますが,必要に応じてインクルードしましょう. 以下にインクルードの一例を載せておきます.

#include <Eigen/Core> /* Coreモジュールをインクルード */
#include <Eigen/Geometry> /* Geometryモジュールをインクルード */
#include <Eigen/Householder> /* Householderモジュールをインクルード */

ご覧のように各モジュールはEigenフォルダ下に収められています. Macな方はBrewでEigenをインストールされることと思われますが,その場合はeigen3フォルダ中にEigenフォルダが格納される可能性がありますので,以下のようなインクルード文が良いでしょう.

#include <eigen3/Eigen/Dense> /* Denseモジュールをインクルード */

もしくは,前者のインクルード文を採用して,GCCの-Iオプションでインクルード先パスを指定する方法もあります.

g++ input.cpp -o output -I /usr/local/include/eigen3/

こうすることでeigen3フォルダ下に必要なライブラリがあることをGCCへ教えることが出来ます. さて,Eigenの基礎知識とコンパイル方法が分かったところで,基本的な使い方を紹介したいと思います.

ベクトル,行列の書き方

Eigen3では,ベクトルの宣言は以下の様に書けます.

Vector2i a;

これでaという名前の2次元ベクトルを宣言できました. Vector2iという部分がクラス名ですが,この"2"が2次元ベクトルであることを意味しています. 3次元ベクトルを宣言したければこの部分を"3"に変えれば良いわけです. また最後の"i"はこのベクトルの成分がInteger,すなわち整数型であることを意味しています. よって,成分がFloat型な3次元ベクトルを宣言したければ以下の様に書けるわけです.

Vector3f b;

これでbという名の3次元ベクトルを宣言しました. ベクトルの成分として扱える型はiとfとd,すなわちInteger型とFloat型とDouble型の3つです. 必要に応じ変更しましょう. さて,次は宣言したベクトルに成分を代入する方法です. 例えば次の3次元ベクトルを定義したいとします. $$ c = \left( \begin{array}{c} 2\\ -1\\ 3 \end{array} \right) $$ cという名のベクトルは,3次元ベクトルで成分はInteger型なので,Vector3iクラスにすれば良さそうです. 1つ目に紹介する方法は,宣言と成分の代入を別々に行う方法です.

Vector3i c;
c << 2, -1, 3;

cを3次元Integer型ベクトルとして宣言してから,"<<"演算子で成分を代入する方法です. これから紹介する行列もこれと同様の記法で成分を代入できます. 2つ目は宣言と同時に代入する方法です.

Vector3i c(2, -1, 3);

ベクトルの定義ならば,これが一番シンプルな書き方です. そして3つ目は成分を番号指定して代入する方法です.

Vector3i c;
c(0) = 2;
c(1) = -1;
c(2) = 3;

配列の様に代入したい成分の位置を番号で指定して代入する方法です. 数学ではベクトルの成分は1番目から始まりますが,Eigen3では配列と同様0番目から数え始めます.

行列は以下の様に宣言します.

Matrix3i A;

これで3行3列の正方行列(3次正方行列)を宣言できました. 行と列の数を変えたい場合は,次の様に宣言します.

MatrixXd B(2, 5);

クラス名において次数を指定していた部分を"X"とし,行列名の横に丸括弧で行数と列数を書きます. これでBという2行5列の行列が宣言されました. なお,この成分にはDouble型を持てます. クラス名の最後がdだからです. 行列成分の代入方法は,ベクトルの代入方法1つ目と同じです. 例として,次の行列Cを定義してみましょう. $$ C = \left( \begin{array}{cc} 1 & 5\\ -2 & 3\\ 4 & -1 \end{array} \right) $$ 行列Cは3行2列であり,その成分は整数ですのでMatrixXiクラスで宣言すれば良さそうです.

MatrixXi C(3, 2);
C << 1,  5,
    -2,  3,
     4, -1;

"<<"演算子を使って行列へ成分を代入して行きます. 見やすさのため改行していますが,改行の必要はありません. 成分はすべてカンマ","で区切ります. さて,ベクトルと行列の定義方法が分かったので,次に簡単な計算を行わせてみましょう.

基本的な計算

ベクトルと行列の基本的な計算を,エクササイズを通して行ってみましょう.

エクササイズ1

2つのベクトルが$ a = \left( \begin{array}{c} 2\\ -1 \end{array} \right),\ b = \left( \begin{array}{c} 3\\ 1 \end{array} \right) $である時,$ 3a-2b $はいくつになるでしょうか?

exercise01.cpp

#include <iostream>
#include <Eigen/Dense>

using namespace std;
using namespace Eigen;

int main()
{
  Vector2i a;
  a << 2, -1;

  Vector2i b;
  b << 3, 1;

  Vector2i solution = 3*a - 2*b;

  cout << "3a - 2b =" << endl << solution << endl;
}

ベクトル同士の積は"*"演算子で,差は"-"演算子で行えます. また,和は"+"演算子で行えます. "ベクトル/スカラ"と書くことでEigenクラスのベクトルをあるスカラ値で割ることが出来ます. スカラはC++の数値リテラル,もしくは変数です.

エクササイズ2

2つのベクトルが$ a = \left( \begin{array}{c} 3\\ -2 \end{array} \right),\ b = \left( \begin{array}{c} 4\\ 3 \end{array} \right) $である時,$ a $と$ b $の内積はいくつでしょうか?

exercise02.cpp

#include <iostream>
#include <Eigen/Dense>

using namespace std;
using namespace Eigen;

int main()
{
  Vector2i a;
  a << 3, -2;
  
  Vector2i b;
  b << 4, 3;

  int solution = a.dot(b);

  cout << "a \u2022 b = " << solution << endl;
}

内積を求めたいベクトルのメソッドdotに被内積ベクトルを与えることで,2つのベクトルの内積を求めることが出来ます. 因みに,エスケープ文字である"\u2022"はドット"・"です.

エクササイズ3

2つのベクトルが$ a = \left( \begin{array}{c} 1\\ 2\\ 3 \end{array} \right),\ b = \left( \begin{array}{c} 2\\ -3\\ 1 \end{array} \right) $である時,$ a $と$ b $の外積はいくつでしょうか?

exercise03.cpp

#include <iostream>
#include <Eigen/Dense>

using namespace std;
using namespace Eigen;

int main()
{
  Vector3i a;
  a << 1, 2, 3;

  Vector3i b;
  b << 2, -3, 1;

  Vector3i solution = a.cross(b);

  cout << "a \u00D7 b = " << endl << solution << endl;
}

外積を求めたいベクトルのメソッドcrossに被外積ベクトルを与えることで,2つのベクトルの外積を求めることが出来ます. 内積とは違い,Geometryモジュールのインクルードが必要になります. エスケープ文字の"\u00D7"はクロス"×"です.

エクササイズ4

2つの行列が$ A = \left( \begin{array}{cc} 3 & 1\\ 8 & 3\\ 1 & 0 \end{array} \right),\ B = \left( \begin{array}{cc} 1 & 5\\ -2 & 3\\ 4 & -1 \end{array} \right) $である時,$ 2A-3B $はいくつになるでしょうか?

exercise04.cpp

#include <iostream>
#include <Eigen/Dense>

using namespace std;
using namespace Eigen;

int main()
{
  MatrixXi A(3, 2);
  A <<
    3, 1,
    8, 3,
    1, 0;

  MatrixXi B(3, 2);
  B <<
    1,  5,
   -2,  3,
    4, -1;

  MatrixXi solution = 2*A - 3*B;

  cout << "2A-3B = " << endl << solution << endl;
}

以上,エクササイズも含めたEigenの紹介でした. 行列をC++で扱うことの出来る素晴らしいライブラリです. ここで紹介したことはほんの一部です. 他にも素晴らしい機能がたくさんあるので,ドキュメンテーションやブログ,フォーラムなどを参考に,様々なことを試してみましょう.

2574973320718972860 https://www.storange.jp/2014/12/eigen.html https://www.storange.jp/2014/12/eigen.html Eigenで線形代数 2014-12-11T15:36:00+09:00 https://www.storange.jp/2014/12/eigen.html Hideyuki Tabata Hideyuki Tabata

200 200

72 72