添付のlibSPARSE.hppを自分のターゲットに追加してください．では，遊んでみましょう．
こいつは，疎行列を勉強しようとしたら＠＃！？格納方式とか面倒だったので，「いやほれ，std::mapって疎行列じゃん．これ，その道のプロが気合いで書いたやつだから，最高性能のはず」ということで，std::mapで計算できるのではないか？と試したものです．オモチャだと思ってほっといたら，平塚くんが装置の温度場解くのに使ってて笑えたやつ
SPARSE::Vector クラス
SPARSE::Vectorは, 疎ベクトルを定義します.
#include "libSPARSE.hpp" ... SPARSE::Vector vec; vec[0]=2.0; //SPARSE::Vectorでは, こうやって定義したところだけ，メモリーを確保します. vec[4]=1.0; vec.Print(6);//6成分ベクトルだと考えて表示してみます.
出力はこんな感じ：
[0] 2.000000[1] -[2] -[3] -[4] 1.000000[5] -
SPARSE::Vector.Print()は, 成分が存在しない場合, - と表示します．
SPARSE::Equationクラス
SPARSE::Equationは, 疎方程式を定義します. 対角成分が0ではないと仮定しているので, 行番号(0,1,...)と対角成分を指定して初期化します．
#include "libSPARSE.hpp" ... SPARSE::Equation eq(3,3.0); //これはi=3の方程式で, 対角成分は3 eq[5]=3.2; //eq[3]=3 以外に0ではない係数を設定できます eq.Print(6);//6成分ベクトルだと考えて表示してみます. eq=0; //スカラー値を代入すると，現状非ゼロの成分だけに代入されます．
出力はこんな感じ：
[0] -[1] -[2] -[3] 3.000000[4] -[5] 3.200000
SPARSE::Vector.Print()は, 成分が存在しない場合, -と表示します．
SPARSE::Coefクラス
SPARSE::Coefは, 疎行列を定義します. 対角成分が0ではないと仮定しているので, 各行に行番号(0,1,...)と対角成分を与えて初期化します.
#include "libSPARSE.hpp" ... SPARSE::Coef coef; for(int i=0;i<6;i++) coef.New(i,1.0);//6x6疎行列を対角成分1で作成 //必要な要素を設定します. もちろん対角要素を再設定しても良いです. coef[0][2]=-1; coef[0][1]=2.3; coef[0][2]=-1; coef[0][3]=+1; coef[0][4]=-3; coef[0][5]=+1; coef[1][0]=-2.3; coef[1][2]=3.0; coef[2][1]=-3.0; coef[2][3]=3.0; coef[3][2]=-3.0; coef[3][4]=3.0; coef[4][3]=-3.0; coef[4][5]=3.0; coef[5][4]=-3.0; coef[5][0]=-1.0; eq.Print(6);//教育用ですので画面に表示することも可能. //でかい行列を表示すると何が起こるかわかんねえだ
出力はこんな感じ：
6x6 ( 0) [0] 1.000000[1] 2.300000[2] -1.000000[3] 1.000000[4] -3.000000[5] 1.000000 ( 1) [0] -2.300000[1] 1.000000[2] 3.000000[3] -[4] -[5] - ( 2) [0] -[1] -3.000000[2] 1.000000[3] 3.000000[4] -[5] - ( 3) [0] -[1] -[2] -3.000000[3] 1.000000[4] 3.000000[5] - ( 4) [0] -[1] -[2] -[3] -3.000000[4] 1.000000[5] 3.000000 ( 5) [0] -1.000000[1] -[2] -[3] -[4] -3.000000[5] 1.000000
SPARSE::Coef.Print()は, 成分が存在しない場合, -と表示します．
SPARSE::CoefとSPARSE::Vectorの掛け算はできます：
auto result=coef*vec;
掛け算を行うと, vecもresultも密ベクトルに化けますので注意してください.
SPARSE::Coefには, 後から列を追加できます:
for(int i=0;i<coef.Dimension();i++) coef[i][coef.Dimension()]=result[i]; coef.Print();
これはつまり\[\text{A}\boldsymbol{x}=\boldsymbol{b}\]を解くために, Aの隣にbを書きたして\[\text{A} \boldsymbol{b}\]という行列にして, こいつを基本変形するって作戦ですよ．ちなみに書き足すbベクトルが複数であると, 複数の連立方程式が解けるだけです．こんな感じの行列になります(この例ではbを2個入れてるな)：
6x8 ( 0) [0] 1.000000[1] 2.300000[2] -1.000000[3] 1.000000[4] -3.000000[5] 1.000000[6] -1.000000[7] 0.000000 ( 1) [0] -2.300000[1] 1.000000[2] 3.000000[3] -[4] -[5] -[6] -4.600000[7] 1.000000 ( 2) [0] -[1] -3.000000[2] 1.000000[3] 3.000000[4] -[5] -[6] 0.000000[7] 4.000000 ( 3) [0] -[1] -[2] -3.000000[3] 1.000000[4] 3.000000[5] -[6] 3.000000[7] 9.000000 ( 4) [0] -[1] -[2] -[3] -3.000000[4] 1.000000[5] 3.000000[6] 1.000000[7] 16.000000 ( 5) [0] -1.000000[1] -[2] -[3] -[4] -3.000000[5] 1.000000[6] -5.000000[7] 25.000000
SPARSE::Coefをコピーすると，非ゼロ要素だけがコピーされます：
auto newMAT=oldMAT;
SPARSE::Coef は, 目次を作成するとLU分解ができます.
coef.MakeIndex(); //目録作成 coef.EraseUp(); //上三角部分を吹っ飛ばす coef.EraseDown(); //下三角部分を吹っ飛ばす coef.ClearIndex(); //目録を消去・・・しとかないと後々まずいので消しといてね coef.Print();
これで対角行列になってしまいました:
6x8 ( 0) [0] 1.000000[1] -[2] -[3] -[4] -[5] -[6] 2.000000[7]-10.424262 ( 1) [0] -[1] 1.000000[2] -[3] -[4] -[5] -[6] -0.000000[7] -4.191884 ( 2) [0] -[1] -[2] 1.000000[3] -[4] -[5] -[6] -0.000000[7] -6.261306 ( 3) [0] -[1] -[2] -[3] 1.000000[4] -[5] -[6] -0.000000[7] -0.771449 ( 4) [0] -[1] -[2] -[3] -[4] 1.000000[5] -[6] 1.000000[7] -3.004156 ( 5) [0] -[1] -[2] -[3] -[4] -[5] 1.000000[6] -0.000000[7] 5.563270
ってこた, ケツの方に解が並ぶので, 取り出して答えだということにします:
SPARSE::Vector solution; for(int i=0;i<coef.Dimension();i++) solution[i]=coef[i][coef.Dimension()]; solution.Print(coef.Dimension());
答えは，あってますね！
[0] 2.000000[1] -0.000000[2] -0.000000[3] -0.000000[4] 1.000000[5] -0.000000
| 添付 | サイズ |
|---|---|
| libSPARSE.hpp__0.zip | 4.78 KB |