「C++/Boost 基礎文法最速マスター」 というサイトのプログラムを試させてもらう。
testboost1.cpp |
#include <iostream> // std::cout を使うのに必要 #include <boost/format.hpp> // boost::format を使うのに必要 int main() { std::cout << boost::format("%s\n") % "Hello, World"; return 0; } |
こんなふうにコンパイル&実行 |
c++ -I/opt/local/include testboost1.cpp ./a.out |
#include <iostream> // std::cout を使うのに必要 #include <boost/format.hpp> // boost::format を使うのに必要 int main() { std::cout << boost::format("%s, %s\n") % "Hello" % "World"; return 0; } |
#include <iostream> // std::cout を使うのに必要 #include <boost/format.hpp> // boost::format を使うのに必要 int main() { std::cout << boost::format("%1%, %2%\n") % "Hello" % "World"; return 0; } |
さて、Boost には、 uBLAS というライブラリィがある。 機能的には BLAS ?
でもちゃんと効率をあげることを考えているのか謎である (C で書いた場合と比較しているページとかあったけれど、 その内容が????????だった)。 普通 BLAS 使おうという人はスピード狂だと思うけど、 そういう人達を納得させる説明が見つからなかった。
以下のプログラムは 「boost::numeric::ublas 線形代数ライブラリの使い方」 にあったものをいじった (元々は A[3][3] のように書いてあったが、 それだとコンパイル出来ないので、A(3,3) のように直す)。
uBLAS にある lu_factorize(), lu_substitute() を試している (なるほど、なるほど)。
testublas.cc |
/* * testublas.cc * g++ -I/opt/local/include testublas.cc */ #include <boost/numeric/ublas/matrix.hpp> // 普通の行列用のヘッダ #include <boost/numeric/ublas/triangular.hpp> // 三角行列用のヘッダ.前進消去,後退代入に必要 #include <boost/numeric/ublas/lu.hpp> // LU分解,前進消去,後退代入用のヘッダ #include <boost/numeric/ublas/io.hpp> // ストリーム入出力用ヘッダ using namespace boost::numeric::ublas; // boost::numeric::ublas 名前空間を使用 #include <iostream> using namespace std; int main() { typedef boost::numeric::ublas::vector<double> dvector; typedef boost::numeric::ublas::matrix<double> dmatrix; // 逆行列を求めるサンプル // 右辺の定数行列を3次単位行列に dmatrix B = identity_matrix<double>(3); // 左辺の係数行列 dmatrix A(3,3); A(0,0) = 1.0; A(0,1) = 2.0; A(0,2) = 3.0; A(1,0) = 2.0; A(1,1) = 1.0; A(1,2) = -3.0; A(2,0) = 4.0; A(2,1) = -3.0; A(2,2) = 2.0; dmatrix CopyOfA(A); //元の行列のコピーを取っておく // LU分解の際の枢軸選択によって起こる行交換の情報が入る行列 permutation_matrix<> pm(A.size1()); // 係数行列と定数ベクトルの表示 cout << "A = " << A << endl; cout << "B = " << B << endl; // LU 分解します. lu_factorize(A,pm); // 行交換を考慮した前進消去と後退代入 lu_substitute(A,pm,B); // 解の表示 cout.precision(4); //小数点以下4桁まで表示 cout << "A の逆行列 X = " << B << endl; cout << "A * X = " << prod(CopyOfA,B) << endl; return 0; } |
コンパイル&実行 |
$ g++ -I/opt/local/include testublas.cc $ ./a.out A = [3,3]((1,2,3),(2,1,-3),(4,-3,2)) B = [3,3]((1,0,0),(0,1,0),(0,0,1)) A の逆行列 X = [3,3]((0.1014,0.1884,0.1304),(0.2319,0.1449,-0.1304),(0.1449,-0.1594,0.04348)) A * X = [3,3]((1,5.551e-17,0),(-1.11e-16,1,8.327e-17),(0,0,1)) $ |
そもそも何で逆行列なんだろう (ネットには逆行列を計算するプログラムがいくつかあった) …まあ、いいや。
桂田 祐史