Motivation
なぜrust?
- 最近C++がつらい
- rustが流行っているらしい
- https://imoz.jp/note/rust-functions.html
環境
- ubuntu 16.04 64bit
rustのインストール
このコマンドで簡単にインストールする
curl https://sh.rustup.rs -sSf | sh
以下のようなメッセージができる． 初心者なのでdefaultのstableをinsatllする．
default host triple: x86_64-unknown-linux-gnu
default toolchain: stable
modify PATH variable: yes
1) Proceed with installation (default)
2) Customize installation
3) Cancel installation
インストールが終わると以下のメッセージが出てくるので
To configure your current shell run source $HOME/.cargo/env
.bashrcかなんかに
source $HOME/.cargo/env
を追加する
rustc --version
を入力して，それらしいものがでてきたらとりあえずOK．
References
crateの作成
Rust では create (木箱) 単位でプログラムを開発するらしい
今回は簡単な行列計算をやってみるので，mat_calcという名前でcrateを作る
cargo new mat_calc --bin
なんかできた
output: Created binary (application) `mat_calc` project
ディレクトリを覗いてみると
$ ls mat_calc . .. Cargo.toml .git .gitignore src
どうやら.gitも生成されているらしいぱっと見た感じ.gitignoreもいい感じに定義されてた
実行テスト
src/main.rsにはすでにHello Worldを出力するプログラムが入っているのでとりあえず実行してみる
とりあえずディレクトリに移動して
$ cd mat_calc $ cargo run Compiling mat_calc v0.1.0 (file:///home/kurenaif/Desktop/mat_calc) Finished dev [unoptimized + debuginfo] target(s) in 0.53 secs Running `target/debug/mat_calc` Hello, world!
無事出力に成功した
References
dependenciesの設定
行列を扱いたいけどFortranでもないしrustにはそんなライブラリは標準搭載されていない. rustで行列が扱えるライブラリを探したら以下の3つが見つかった
- https://github.com/sebcrozet/nalgebra
- https://github.com/termoshtt/ndarray-linalg
- https://github.com/masonium/linxal
一番使いやすそうで，commitが盛んなnalgebraを使うことにした．
速度面はよくわからなかった．
ここで, Cargo.tomlを編集し，nalgebraを使用できるようにする(pip installに近いイメージ)
[dependencies] nalgebra = "0.11.0"
versionは2017-05-05の時点での最新バージョンを使用した.
githubなりなんなりをみてうまいこと指定する．
githubのURLでの指定も可能だけど https://crates.io に登録されているのでこんな感じでURLなしでもcrateを入れることが出来る.
References
連立方程式を解いてみる
行列の表示
以下のソースコードをsrc/main.rsに実装し，実行してみる
今後の流体計算を考慮して普通のMatrixよりも実行時に長さを決められるDMatrixのほうが都合がよさそうなのでDMatrixを使用する.
使い方その他もろもろはソースコードみて察していただきたい．
個々で気をつけなければならないのが行列の数字はすべて実数で表記しないとちゃんとprintlnしてくれない.
つまり 1,2,3… のような整数を書く場合でも, 1., 2., 3. と書いたほうが後々楽になる．
rustは整数と小数，また同じ整数どうしや小数どうしてもbit幅が違えばコンパイルエラーを出すほど厳密．
extern crate nalgebra as na; fn main(){ let mat = na::DMatrix::from_iterator(3,3,[ 1.,2.,3., 4.,5.,6., 7.,8.,9. ].iter().cloned()); println!("{}",mat); let b = na::DVector::from_iterator(3,[ 1.,2.,3. ].iter().cloned()); println!("{}",b); }
$ cargo run ┌ ┐ │ 1.000 4.000 7.000 │ │ 2.000 5.000 8.000 │ │ 3.000 6.000 9.000 │ └ ┘ ┌ ┐ │ 1.000 │ │ 2.000 │ │ 3.000 │ └ ┘
- めっちゃいい感じの出力をしてくれる.(DMatrixじゃなくて普通のMatrixはこうはいかない)
- どうやらソースコードに書いたものに比べて転置したものが出てくるらしい 縦ベクトルが3つ横に並んでるイメージ.
逆行列
方程式を解くためには逆行列が必要不可欠 さっき用意したmatは実は逆行列が存在しない．
面白そうなのでさっきのmatの逆行列を計算してみる．
nalgebraでは，.try_inverse()を用いて逆行列を計算する.
無理だったら無理なりの挙動をしてくれる.
try_inverseを出力するプログラム．
println!("{}") と println!("{:?}")では呼びされるメソッドが微妙に異なっており，後者の方がちょっと詳しいメッセージが出てくることが多い(デバッグ用出力?)．
今回，前者を使うとエラーを吐くので後者を使っている．
extern crate nalgebra as na; fn main(){ let mat = na::DMatrix::from_iterator(3,3,[ 1.,2.,3., 4.,5.,6., 7.,8.,9. ].iter().cloned()); println!("{:?}",mat.try_inverse()); }
$ cargo run None
Noneが出力された．
リファレンスに書いているが，try_inverse()をしてなかった時はNoneが返ってくるらしい
では逆行列が存在する時はどんな出力をするのだろうか?
今回はtry_inverseの結果を"{}"と"{:?}"の二つを出力している．
最初のいもすさんの記事でも書いてある通り，rustでは基本所有権がmoveしちゃうので最初のprintln!でmatの所有権が切れる．
そこで，clone()することで所有権を継続させている．
extern crate nalgebra as na; fn main(){ let mat = na::DMatrix::from_iterator(3,3,[ 1.,3.,5., 2.,3.,8., 2.,1.,5. ].iter().cloned()); println!("{:?}",mat.clone().try_inverse()); println!("{}",mat.try_inverse().unwrap()); }
Some(Matrix { data: MatrixVec { data: [1.4, -2, 1.8, 1.2, -1, 0.4, -0.8, 1, -0.6], nrows: Dynamic { value: 3 }, ncols: Dynamic { value: 3 } }, _phantoms: PhantomData })
┌ ┐
│ 1.400 1.200 -0.800 │
│ -2.000 -1.000 1.000 │
│ 1.800 0.400 -0.600 │
└ ┘
一つ目の出力がtry_inverse()の中身を詳細に出力したものである．
Some()がついていて，実はこのままでは{}の出力ができない．
ここでSome()は値がNullになるかもしれない的なやつらしい．
Some()を取るためには，unwrap()をすると良い．
詳しくはここを参照する:
これの結果があってるかわからないので元のmatにかけて単位行列が出ること確認する
rustでは演算子のオーバーロードができる．
このライブラリでは行列の掛け算もちゃんと定義されているので，以下のソースコードみたいな感じで表記できる.
extern crate nalgebra as na; fn main(){ let mat = na::DMatrix::from_iterator(3,3,[ 1.,3.,5., 2.,3.,8., 2.,1.,5. ].iter().cloned()); println!("mat\n{}",mat.clone()); println!("inv(mat)\n{}",mat.clone().try_inverse().unwrap()); println!("mat*inv(mat)\n{}",mat.clone().try_inverse().unwrap()*mat); }
$ cargo run mat ┌ ┐ │ 1.000 2.000 2.000 │ │ 3.000 3.000 1.000 │ │ 5.000 8.000 5.000 │ └ ┘ inv(mat) ┌ ┐ │ 1.400 1.200 -0.800 │ │ -2.000 -1.000 1.000 │ │ 1.800 0.400 -0.600 │ └ ┘ mat*inv(mat) ┌ ┐ │ 1.000 -0.000 0.000 │ │ 0.000 1.000 0.000 │ │ 0.000 0.000 1.000 │ └ ┘
ちゃんと逆行列が計算されているっぽい
連立方程式を解いてみる
いよいよ本題
とりあつかう問題はこれの例1: 連立方程式の解き方
連立方程式を $$ x = A^{-1} b$$ の形にしてinvして解くだけ
extern crate nalgebra as na; fn main(){ let a = na::DMatrix::from_iterator(3,3,[ 1.,1.,0., 2.,-1.,1., 3.,2.,1. ].iter().cloned()); println!("A\n{}",a); let b = na::DVector::from_iterator(3,[ 0.,3.,-1. ].iter().cloned()); println!("B\n{}",b); let invA = a.try_inverse().unwrap(); println!("invA\n{}",invA); println!("x\n{}", invA*b) }
A ┌ ┐ │ 1.000 2.000 3.000 │ │ 1.000 -1.000 2.000 │ │ 0.000 1.000 1.000 │ └ ┘ B ┌ ┐ │ 0.000 │ │ 3.000 │ │ -1.000 │ └ ┘ invA ┌ ┐ │ 1.500 -0.500 -3.500 │ │ 0.500 -0.500 -0.500 │ │ -0.500 0.500 1.500 │ └ ┘ x ┌ ┐ │ 2.000 │ │ -1.000 │ │ 0.000 │ └ ┘
ページの答えと一致してるしちゃんと解けてるっぽい．
おまけ
今回色々やったリポジトリをここにおいておきます github.com
次回は簡単な熱拡散問題とか解いてみるかも
