Eigen ー C++で線形代数を!(2)

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

(このエントリーは旧ブログからの転載です)

ここではC++用線形代数ライブラリ eigen を個人的に習得するために作ったサンプルプログラムを公開しています(一部、本家のチュートリアルと重複しています)。使ってみて感じたのは C++ という言語で ruby や R 的な配列の柔軟さを実現していてすごいなーということでした。ブロードキャストの例に出てくるメソッドチェーンなんて結構感動(式テンプレートのおかげで、おそらくパフォーマンスは落ちない)。すばらしいライブラリです。以下、ライブラリの基本的な説明です。固有値計算などいわゆる線形代数特有の演算についてはいずれ書きます。

【追記】 続きを書きました。Eigen ー C++で線形代数を!(3)

目次

注意!

以下のすべてのサンプルでは表記の簡略化のために、以下のコードを省略しています。Eigen ライブラリの読み込み & 行列表示と、関数名表示のための簡単なマクロ & Eigen 名前空間。

#include <Eigen/Core>
#define PRINT_MAT(X) cout << #X << ":\n" << X << endl << endl
#define PRINT_MAT2(X,DESC) cout << DESC << ":\n" << X << endl << endl
#define PRINT_FNC    cout << "[" << __func__ << "]" << endl

using namespace Eigen;

使用バージョンは Eigen 3.0 beta-2 です。Eigen はヘッダのみで動作するライブラリなので、オブジェクトファイルへのリンクなどは必要ありません。

基本の操作(1) - 行列・ベクトルの初期化、コピー、サイズ変更

固定サイズ行列・動的なサイズの行列
固定サイズベクトル・動的サイズベクトル
行列の初期化 (1) : カンマオペレータ
行列の初期化 (2) : 特殊な行列形式
行列の初期化 (3) : 要素の直接参照
ベクトルの初期化
行列のメモリ配置の変更
行列のコピーとスワップ
行の入れ替え、列の入れ替え
行列のリサイズ

基本の操作(2) - 要素・行・列・ブロック行列(部分小行列)への参照

行列要素へのアクセス
行列の行・列へのアクセス
行列の部分行列へのアクセス (1) : 部分行列の取り出し
行列の部分行列へのアクセス (2) : 部分行列への代入

基本の操作(3) - ベクトル・行列情報の取得

ベクトルと行列のサイズの取得

基本の操作(4) - 算術演算

行列の要素ごとの演算
行列演算 (1) : 行列同士の和、差、積
行列演算 (2) : 行列とスカラーの和、差、積

基本の操作(5) - 比較演算

比較演算 (1) : 行列の各要素に対して条件判定 → bool型行列を返す
比較演算 (2) : 行列の各要素に対して条件判定 → 結果を集約して返す

基本の操作(6) - 対角行列・転置行列

対角行列
転置行列

基本の操作(7) - 情報集約、ノルム、ビジター

行列の情報集約:要素ごとの和、積、平均、最大、最小、トレース
行列の部分的情報集約:行ごと・列ごとの和、積、平均、最大、最小
行列ノルム、ベクトルノルム
ビジター:最大、最小要素の添え字を得る

基本の操作(8) - ブロードキャスティング

ブロードキャスト (1) :行、列に対してある操作を一気に行う
ブロードキャスト (2) :行、列に対してある操作を一気に行う

基本の操作(1) - 行列・ベクトルの初期化、コピー、サイズ変更

固定サイズ行列・動的なサイズの行列

eigen ではさまざまな名前の行列クラス、ベクトルクラスが存在しますが、命名規則を覚えてしまえば簡単です。

たとえば、CGなどで使いそうなdouble型の三次元行列クラスは「Matrix3d」です。「3」が次元をあらわしています。他方「Matrix3d」の「d」は「dimension」ではなく、「double」です。float なら Matrix3f だし、int なら Matrix3i となります。bool は b、complex は c です。

一方、次元のほうは4次元まで指定できます。それ以上の場合はテンプレート引数で指定します。テンプレートは Matrix<double, 6, 6> のような感じです(double型6行6列の場合)。さらにコンパイル時にはサイズが不明な場合のために、動的な行列も用意されています。それらは MatrixXd もしくは Matrix<double, Dynamic, Dynamic> で宣言します。

固定サイズベクトル・動的サイズベクトル

ベクトルに関しても行列同様、多くの型がありますが、行列とほとんど同じ命名規則です。MatrixXx が VectorXx になるだけですね。

行列の初期化 (1) : カンマオペレータ

Eigen ではさまざまな行列の初期化の方法があります。このカンマオペレータはほかのライブラリでは見かけませんが、非常に楽な初期化の仕方です。

void Matrix_Init_Test1()
{
    /* 行列クラスの初期化 */
    Matrix3i A;  // 3x3 の int 型行列

    // (1) カンマ演算による初期化
    A << 1,2,3,
         4,5,6,
         7,8,9;
 
    PRINT_MAT(A);
}

行列の初期化 (2) : 特殊な行列形式

行列を特殊行列で初期化します。このとき注意すべきは右辺と左辺で同じ型を指定しなければならない、ということです。

void Matrix_Init_Test2(int n){
    // (2) 特殊型のベクトルで割り当て
    int n = 4;
    MatrixXd A1 = MatrixXd::Zero(n,n);          // すべて0
    MatrixXd A2 = MatrixXd::Ones(n,n);          // すべて1
    MatrixXd A3 = MatrixXd::Constant(n,n,2);    // 定数ベクトル
    MatrixXd A4 = MatrixXd::Random(n,n);        // ランダム
    MatrixXd A5 = MatrixXd::Identity(5,5) // 単位行列(正方行列でなくてもOK)

    Matrix3d B = Matrix3d::Identity()  // 固定サイズのときはサイズ指定必要なし。

    PRINT_MAT(A1);
    PRINT_MAT(A2);
    PRINT_MAT(A3);
    PRINT_MAT(A4);
    PRINT_MAT(A5);
}

行列の初期化 (3) : 要素の直接参照

要素への参照を添え字を通して行う方法です。直感的にわかりやすい方法です。上の二つの方法で済む場合は、バグをなくすためにも使わないほうがよいのかもしれません。

void Matrix_Init_Test3(){
    MatrixXd A(3,3);
    for(int i=0;i<3;++i){
        for(int j=0;j<3;++j) A(i,j)=i+j;
    }
    PRINT_MAT(A);
}

ベクトルの初期化

ベクトルの初期化に関しても行列の場合と同じだと思ってよいでしょう。

void Vector_Init_Test()
{
    /* ベクトルクラスの初期化
     * 行列クラスと同様に配列サイズの固定バージョン Vector3d などがある */
    VectorXd u(4);

    // (1) カンマ演算による初期化
    u << 1,2,3,4;
 
    PRINT_MAT(u);
 
    // (2) 特殊型のベクトルで割り当て
    int n = 4;
    VectorXd v1 = VectorXd::Zero(n);          // すべて0
    VectorXd v2 = VectorXd::Ones(n);          // すべて1
    VectorXd v3 = VectorXd::Constant(n,2);    // 定数ベクトル
    VectorXd v4 = VectorXd::Random(n);        // ランダム
    VectorXd v5 = VectorXd::LinSpaced(n,1,2); // 線形補間
    VectorXd v6 = VectorXd::Unit(n,2);        // 単位ベクトル
 
    PRINT_MAT(v1);
    PRINT_MAT(v2);
    PRINT_MAT(v3);
    PRINT_MAT(v4);
    PRINT_MAT(v5);
    PRINT_MAT(v6);
}

行列のメモリ配置の変更

行列のメモリ配置を行優先にするか、列優先にするかです。大規模でシビアな計算では必要になってきます。省略可能なテンプレート引数を明示的に指定することで変更できますが、通常はデフォルトでよいでしょう。あとは下のサンプルのなかのコメントを見てください。

// Row-major/Column-major test
void Matrix_RowCol_Major_Test()
{
    PRINT_FNC;

    Matrix<double,2,3,rowmajor>; Arow;
    Matrix<double,2,3,colmajor> Acol;
    Arow <<
        1,2,3,
        4,5,6;
    Acol <<
        1,2,3,
        4,5,6;

    // 行列の表示
    PRINT_MAT(Arow);
    PRINT_MAT(Acol);
 
    // ストレージの格納順に表示
    cout << "Storage ordering of Arow: " << endl;
    for(int i=0;i<6;++i){
        cout << Arow(i) << " ";
    }
    cout << endl;
 
    cout << "Storage ordering of Acol: " << endl;
    for(int i=0;i<6;++i){
        cout << Acol(i) << " ";
    }
    cout << endl;
    /* 注意!
     * デフォルトは ColMajor になっている.
     * カンマオペレータによる初期化は「ストレージの方向」に依存しない!
     * (常に横方向に伸びる,見た目が直感的になるから?)
     * ストレージの格納方法が異なることは上の Print storage の部分から
     * 確認できる.
     */
    cout << endl;
}

行列のコピーとスワップ

行列の代入は参照ではなく、要素のコピーが発生します。したがって関数に渡すときは基本的に参照渡しにすべきです。行列オブジェクト同士のスワップは swap メンバー関数を使うことで、要素のコピーを発生させることなく実行可能です。

void Matrix_Copy_Swap_Test()
{
    /* コピーおよびスワップ(交換)*/
    PRINT_FNC;

    MatrixXd A_orig = MatrixXd::Random(3,3);
    MatrixXd A_dump = A_orig;  /* ディープコピー(参照でない!) */

    A_orig(0,0) = 1111;

    /* 一方を変更しても他方には影響しない */
    PRINT_MAT(A_orig);
    PRINT_MAT(A_dump);

    /* スワッピング! */
    cout << "  (Swapping A_orig <--> A_dump)" << endl << endl;
    A_orig.swap(A_dump);
    PRINT_MAT(A_orig);
    PRINT_MAT(A_dump);
}

行の入れ替え、列の入れ替え

行列クラスの col(int), row(int) メンバー関数を使うと列ベクトル、行ベクトルへの参照が得られます。これと swap を組み合わせることで行、列の入れ替えが可能になります。

void Matrix_RowCol_Swap_Test()
{
    MatrixXd A = MatrixXd::Random(3,3);
    PRINT_MAT(A);

    A.col(0).swap(A.col(2));
    PRINT_MAT2(A,"Swapping A[,0] <--> A[,2]");
}

行列のリサイズ

サイズ変更です。コメントにもあるとおり、別のヒープ領域にメモリを確保し、値をコピーしないため、オリジナルの要素の値が(サイズを縮小した場合にさえ)保持されないことに注意してください。

void Matrix_Resize_Test()
{
    PRINT_FNC;

    MatrixXd A = MatrixXd::Random(3,3);
    cout << "[Original]" << endl;
    PRINT_MAT(A);

    A.resize(2,2);
    cout << "[After resize (smaller)]" << endl;
    PRINT_MAT(A);

    A.resize(4,4);
    cout << "[After resize (bigger)]" << endl;
    PRINT_MAT(A);

    /* 注意!
     * resize すると(サイズを縮小した場合でさえ)もとの
     * 要素は保存されない!
     */
}

基本の操作(2) - 要素・行・列・ブロック行列(部分小行列)への参照

行列要素へのアクセス

void Matrix_Element_Access_Test()
{
    /* 行列成分は [] でなく () 演算子でアクセス可能
     * これらの演算子は「配列添字の妥当性チェック」を行う
     * coeff(i,j) や coeffRef(i,j) で「添字チェックなし」の
     * アクセスが可能になる!
     */
    PRINT_FNC;

    MatrixXd A = MatrixXd::Identity(3,3);
    cout << "A(0,0) = " << A(0,0) << endl;;
    cout << "A.coeff(1,2) = " << A.coeff(1,2) << endl;
    cout << "A.coeffRef(1,2) = " << A.coeff(1,2) << endl;

    A.coeffRef(2,2) = 100;  // OK
    //A.coeff(2,2) = 100;  // NG

    cout << endl;
}

行列の行・列へのアクセス

コメントのとおりです。取り出した参照はベクトルのように扱うことができます。参照なので当然代入も可能。

void Matrix_RowCol_Access_Test()
{
    /* .row(), .col() で行,列への参照を取得できる */
    PRINT_FNC;

    MatrixXd A = MatrixXd::Identity(4,4);
    A.col(0) = VectorXd::Random(4);
    A.row(0) = VectorXd::Random(4);

    PRINT_MAT(A);

    VectorXd c2 = A.col(2);
    PRINT_MAT2(c2, "3rd-col of A");
}

行列の部分行列へのアクセス (1) : 部分行列の取り出し

こちらは部分行列への参照を返します。数学の論文とか教科書とかは結構ブロック行列の形式で書かれていることも多いので、これは非常に便利。

void Matrix_Block_Access_Test1()
{
    /* ブロック行列を得る */
    MatrixXd A = MatrixXd::Random(5,5);
    MatrixXd Apart = A.block(0,0,2,2); // (0,0) 成分からはじまる 2x2 小行列
    MatrixXd Apart2 = A.block<2,2>(0,0); // サイズ固定
 
    PRINT_MAT(A);
    PRINT_MAT(Apart);
    PRINT_MAT(Apart2);
}

行列の部分行列へのアクセス (2) : 部分行列への代入

void Matrix_Block_Access_Test2()
{
    /* block(...) で帰ってくる部分行列はもとの行列への
     * 参照になっているので代入による変更ができる
     * 論文などではブロック行列形式で書かれていることも
     * 多いので,これを使いこなすとバグが減りそう.
     */
    MatrixXd A = MatrixXd::Random(5,5);
    A.block(0,0,2,2) = Matrix2d::Identity();

    PRINT_MAT(A);
}

基本の操作(3) - ベクトル・行列情報の取得

ベクトルと行列のサイズの取得

void Vector_Matrix_GetSize_Test(int nrow, int ncol)
{
    /* ベクトル,行列のサイズを得る */
    PRINT_FNC;

    MatrixXd B(nrow, ncol);
    VectorXd v(nrow);

    cout << "Size of vector v is : " << v.size() << endl;
    cout << "Size of matrix B is : ("
         << B.rows() << "," << B.cols() << ")"
         << endl << endl;
}

基本の操作(4) - 算術演算

行列の要素ごとの演算

version 2 と version 3 で API が違います。version 3 では array というメンバ関数を使います。あれこれ説明するより下のサンプルを見たほうが早いと思います。

void Matrix_cwise_func_Test1()
{
    PRINT_FNC;

    MatrixXd A = MatrixXd::Random(3,3);
    PRINT_MAT(A);

    MatrixXd sin_A = A.array().sin();
    PRINT_MAT2(sin_A,"sin(A)");

    MatrixXd cos_A = std::cos(A.array()); // もう一つの方法(より自然?)
    PRINT_MAT2(cos_A,"cos(A)");

    MatrixXd exp_A = A.array().exp();
    PRINT_MAT2(exp_A,"exp(A)");

    MatrixXd log_exp_A = exp_A.array().log();
    PRINT_MAT2(log_exp_A,"log(exp(A))");

    /* 注意!  --- array() は非常に重要!
     * version 2 では A.cwise().sin() などとしていたものが
     * version 3 では A.array().sin() というように変更された.
     * std::sin(A) でも同じこと.
     */
}

行列演算 (1) : 行列同士の和、差、積

非常に簡単にかけます。式テンプレートが使われているということなので、なるべく明示的な一時変数は減らすようにしましょう。

void Matrix_Arith_OP_Test1()
{
    /* 行列同士の +, -, * 演算
     * "*" 演算は行列の意味(coef-wise ではない!) */
    PRINT_FNC;

    MatrixXd A1 = MatrixXd::Random(3,3);
    MatrixXd A2 = MatrixXd::Random(3,3);

    PRINT_MAT(A1);
    PRINT_MAT(A2);

    MatrixXd A3(3,3);

    A3 = A1 + A2;
    PRINT_MAT2(A3,"A1+A2");

    A3 = A1 - A2;
    PRINT_MAT2(A3,"A1-A2");

    A3 = A1 * A2;
    PRINT_MAT2(A3,"A1*A2");

}

行列演算 (2) : 行列とスカラーの和、差、積、商

数学では行列のスカラー倍は定義されるけど、行列とスカラーの和、差は定義されていません。しかし、eigen では A+s*Ones という演算を省メモリで実行するために、行列とスカラーの和、差が定義されています。

void Matrix_Arith_OP_Test2()
{
    /* 行列とスカラーの演算 */
    PRINT_FNC;

    MatrixXd A = MatrixXd::Identity(3,3);

    MatrixXd B = A.array() + 1.0;  // 和
    PRINT_MAT2(B,"A+1.0");

    B = A.array() - 1.0;  // 差
    PRINT_MAT2(B,"A-1.0");

    B = A.array() * 10.0; // 積
    PRINT_MAT2(B,"A*10.0");

    B = A.array() / 10.0; // 商
    PRINT_MAT2(B,"A/10.0");

    // B = 10.0 / A.array();  // <- これはエラー!代わりに次の行のように
    B = A.array().inverse() * 10.0;  // good
    // B = MatrixXd::Constant(A.rows(),A.cols(),10.0).array()/A.array();
    // ↑ bad : たぶんメモリの無駄
    PRINT_MAT2(B,"10.0/A");
}

基本の操作(5) - 比較演算

比較演算 (1) : 行列の各要素に対して条件判定 → bool型行列を返す

void Matrix_Comp_OP_Test1()
{
    // 行列の各成分に対して条件判定して bool 型の行列を返す
    PRINT_FNC;

    MatrixXd A = MatrixXd::Identity(3,3);

    Matrix&lt;bool, dynamic, dynamic&gt; B = A.array()&lt;1.0; // 比較

    PRINT_MAT(A);
    PRINT_MAT2(B, "A<1.0");
}

比較演算 (2) : 行列の各要素に対して条件判定 → 結果を集約して返す

void Matrix_Comp_OP_Test2()
{
    /* 行列の各要素が条件をみたすかどうかチェック
     * (all) すべての要素に対して満たすか?
     * (any) どれかひとつの要素に対して満たすか?
     * (count) 条件を満たす要素数
     */
    PRINT_FNC;

    MatrixXd A = MatrixXd::Identity(3,3);

    PRINT_MAT(A);
    cout << "Does A(i,j)<1.0 hold for all elements? : " << (A.array()<1.0).all() << endl;
    cout << "Does A(i,j)<2.0 hold for all elements? : " << (A.array()<2.0).all() << endl;
    cout << "Does A(i,j)<1.0 hold for any elements? : " << (A.array()<1.0).any() << endl;
    cout << "How many elements do A(i,j)<1.0 hold?  : " << (A.array()<1.0).count() << endl;

    cout << endl;
}

基本の操作(6) - 対角行列・転置行列

対角行列

void Matrix_Diagonal_Init_Test()
{
    /* ベクトルから対角行列を生成する */
    PRINT_FNC;

    VectorXd v = VectorXd::Random(5);
    MatrixXd D = v.asDiagonal();

    PRINT_MAT(v);
    PRINT_MAT2(D,"diag(v)");
}

転置行列

void Matrix_Transpose_Test()
{
    MatrixXd A = MatrixXd::Random(3,4);

    PRINT_MAT(A);

    PRINT_MAT2(A.transpose(),"A^T");
}

基本の操作(7) - 情報集約、ノルム、ビジター

行列の情報集約:要素ごとの和、積、平均、最大、最小、トレース

void Matrix_Reduction_Test()
{
    /* 集約関数を用いて行列の成分を集約する */
    PRINT_FNC;
    MatrixXd A = MatrixXd::Identity(3,3);
    PRINT_MAT(A);

    cout << "Sum  of A(i,j) : " << A.sum() << endl;
    cout << "Prod of A(i,j) : " << A.prod() << endl;
    cout << "Mean of A(i,j) : " << A.mean() << endl;
    cout << "Min  of A(i,j) : " << A.minCoeff() << endl;
    cout << "Max  of A(i,j) : " << A.maxCoeff() << endl;
    cout << "Trace of A(i,j) : " << A.trace() << endl;  // tr(A) は行列の対角成分の和
    cout << endl;
}

行列の部分的情報集約:行ごと・列ごとの和、積、平均、最大、最小

void Matrix_Partial_Reduction_Test()
{
    /* 集約関数を用いて行列の成分を「部分的に=行,列ごとに」集約する */
    PRINT_FNC;

    MatrixXd A = MatrixXd::Identity(3,3);
    PRINT_MAT(A);

    cout << "Sum  of each row of A : \n" << A.rowwise().sum() << endl;
    cout << "Prod of each row of A : \n" << A.rowwise().prod() << endl;
    cout << "Mean of each col of A : \n" << A.colwise().mean() << endl;
    cout << "Min  of each col of A : \n" << A.colwise().minCoeff() << endl;
    cout << "Max  of each row of A : \n" << A.rowwise().maxCoeff() << endl;
    cout << "Notice!!: Result distinguishes between row-major and col-major!!" << endl;
    cout << endl;

    /* 代入も可能 → 結果はすべて縦ベクトル! */
    VectorXd sumAcol = A.colwise().sum();
    VectorXd meanArow = A.rowwise().mean();

    PRINT_MAT2(sumAcol, "Sum of cols of A");
    PRINT_MAT2(meanArow, "Mean of rows of A");
    cout << "Notice!!: Result does not distinguish between row-major and col-major!!" << endl;
}

行列ノルム、ベクトルノルム

void Matrix_Norm_Test()
{
    /* ベクトル・行列のノルムを得る */
    PRINT_FNC;

    VectorXf v(2);
    MatrixXf A(2,2);
    v << -1, 2;
    A << 1,-2,
        -3, 4;

    PRINT_MAT(v);
    cout << "v.squaredNorm() = " << v.squaredNorm() << endl;  // Not rooted!
    cout << "v.norm() = " << v.norm() << endl;
    cout << "v.lpNorm<1>() = " << v.lpNorm<1>() << endl;  // l1-norm (sum of absolute value)
    cout << "v.lpNorm<infinity>() = " << v.lpNorm<infinity>() << endl; // max(abs(coeff))
    cout << endl;

    PRINT_MAT(A);
    cout << "A.squaredNorm() = " << A.squaredNorm() << endl;  // Not rooted!
    cout << "A.norm() = " << A.norm() << endl;
    cout << "A.lpNorm<1>() = " << A.lpNorm<1>() << endl;  // sum of absolute value
    cout << "A.lpNorm<infinity>() = " << A.lpNorm<infinity>() << endl; // max(abs(coeff))
    cout << endl;
}

ビジター:最大、最小要素の添え字を得る

ビジターなるものは私も良くわかっていませんが、チュートリアルにあった、最大、最小要素の添え字を得る、という例はなかなか便利。R でいうところの which.max, which.min ですな。

void Matrix_Visitors_Test(){
    /* どの要素が最小,最大か? (which.min and which.max ?)
     * このようなものを Eigen ではビジターと呼んでいる */
    PRINT_FNC;

    MatrixXd::Index maxRow, maxCol, minRow, minCol;
    MatrixXd A = MatrixXd::Random(4,4);

    //double Amax = A.array().max(&maxRow, &maxCol);  // これはなぜかエラー
    //double Amin = A.array().min(&minRow, &minCol);  // これもなぜかエラー
    double Amax = A.maxCoeff(&maxRow, &maxCol);
    double Amin = A.minCoeff(&minRow, &minCol);

    PRINT_MAT(A);
    cout << "Max of A is " << Amax << " at " << maxRow << ", " << maxCol << endl;
    cout << "Min of A is " << Amin << " at " << minRow << ", " << minCol << endl;
    cout << endl;
}

基本の操作(8) - ブロードキャスティング

ブロードキャスト (1) :行、列に対してある操作を一気に行う

たいそうな名前が付いていてビビリますが、次のようなことらしい。

void Matrix_Broadcast_Test1()
{
    /* 行,列に対してある操作を一気に行う */
    PRINT_FNC;

    MatrixXd A = MatrixXd::Identity(3,3);
    VectorXd v(3);
    v << 1,2,3;

    PRINT_MAT(A);
    PRINT_MAT(v);

    A.colwise() += v;  // すべての列に v を加える <--> rowwise()

    PRINT_MAT2(A,"A+(v,v,v)");
}

ブロードキャスト (2) :行、列に対してある操作を一気に行う

ブロードキャストを使うとこの例のように ruby 的なメソッドチェーンがより柔軟に行えます、という例。

void Matrix_Broadcast_Test2()
{
    /* 行,列に対してある操作を一気に行う */
    PRINT_FNC;

    MatrixXd A = MatrixXd::Identity(3,3);
    VectorXd v(3);
    v << 1,2,3;

    PRINT_MAT(A);
    PRINT_MAT(v);

    /* 行列 A の各列にベクトル v を加えて行方向の最小値をベクトルとして得る */
    PRINT_MAT2((A.colwise() + v).rowwise().minCoeff(),"Rowwise minimum coeff. of (A + (v,v,v))");
}
線形代数と数値解析 (理工系の数学教室) 線形代数と数値解析 (理工系の数学教室)
河村 哲也

朝倉書店
売り上げランキング : 657400

Amazonで詳しく見る

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