【Rust】Rustで数値計算プロジェクトを試す その7

クロージャを渡す実装をしました。 t の引数は結局使ってないので削除しました。

doc.rust-jp.rs

の記法はトレイト境界というもので、これが無いと以下のようなエラーが出るのですが、

error[E0277]: the size for values of type `(dyn std::ops::Fn(f64, f64) -> f64 + 'static)` cannot be known at compilation time
  --> src/main.rs:31:32
   |
31 | fn runge_kutta(x: f64, v: f64, f: Fn(f64, f64) -> f64) -> (f64, f64) {
   |                                ^ doesn't have a size known at compile-time

…なのですが、今のところ意味が分からずに書いてるので後でチェックします。

use gnuplot::{Figure, Caption, Color};


const K: f64 = 1.0;
const M: f64 = 1.0;
const A: f64 = 2.0;
const X0: f64 = 0.0;
const V0: f64 = 1.0;
const DELTA_T: f64 = 0.01;


fn main() {
  let (mut x, mut v) = (X0, V0);
  let mut xs = vec![];
  let mut ts = vec![];

  for i in 1..1000 {
    let t = (i as f64) * DELTA_T;
    let (dx, dv) = runge_kutta(x, v, |x, v| {
      - K * x - A * v
    });
    x += dx;
    v += dv;
    ts.push(t);
    xs.push(x);
  }

  plot(ts, xs);
}

fn runge_kutta<F>(x: f64, v: f64, f: F) -> (f64, f64) 
    where F : Fn(f64, f64) -> f64 {
  let k1v = f(x, v) * DELTA_T / M;
  let k1x = v * DELTA_T;
  let k2x = (v + k1v / 2.0) * DELTA_T; 
  let k2v = f(x + k1x / 2.0, v + k1v / 2.0) * DELTA_T / M;
  let k3v = f(x + k2x / 2.0, v + k2v / 2.0) * DELTA_T / M;
  let k3x = (v + k2v / 2.0) * DELTA_T;
  let k4v = f(x + k3x, v + k3v) * DELTA_T / M;
  let k4x = (v + k3v) * DELTA_T;
  let v = (k1v + 2.0 * k2v + 2.0 * k3v + k4v) / 6.0;
  let x = (k1x + 2.0 * k2x + 2.0 * k3x + k4x) / 6.0;
  (x, v)
}

fn plot(x: Vec<f64>, y: Vec<f64>) {
  let mut fg = Figure::new();
  fg.set_terminal("pngcairo", "output.png");
  fg.axes2d()
    .lines(&x, &y, &[Caption("A line"), Color("black")]);
  fg.show();
}

【Rust】Rustで数値計算プロジェクトを試す その6

ひとまず、Rust習得の第二の中間目標として「Lotka-Volterra方程式を4次のRunge-Kutta法で近似計算する」ことを置いてみます。いきなり2重振り子を計算するのは、極座標とかラグランジアンみたいな概念が必要っぽくて、「大学時代に触れた覚えがあるけど急に思い出せない」ので後にします。

ひとまずRunge-Kutta法の実装です。

shimaphoto03.com

こちらの調和振動子の実装をパクりました。

qiita.com

所有権周りで苦労すると思ってましたが、Copy Traitが実装されているおかげで、関数に渡した後でも使えるので詰まることはありませんでした。

doc.rust-jp.rs

汎用性をもたせるために、関数ポインタとかいろいろ試してみたいですね。

use gnuplot::{Figure, Caption, Color};


const K: f64 = 1.0;
const M: f64 = 1.0;
const X0: f64 = 0.0;
const V0: f64 = 1.0;
const DELTA_T: f64 = 0.01;


fn main() {
  let (mut x, mut v) = (X0, V0);
  let mut xs = vec![];
  let mut ts = vec![];

  for i in 1..1000 {
    let t = (i as f64) * DELTA_T;
    let (dx, dv) = runge_kutta(x, v, t);
    x += dx;
    v += dv;
    ts.push(t);
    xs.push(x);
  }

  plot(ts, xs);
}

fn runge_kutta(x: f64, v: f64, t: f64) -> (f64, f64) {
  let k1v = f(x, v, t) * DELTA_T / M;
  let k1x = v * DELTA_T;
  let k2x = (v + k1v / 2.0) * DELTA_T; 
  let k2v = f(x + k1x / 2.0, v + k1v / 2.0, t + DELTA_T / 2.0) * DELTA_T / M;
  let k3v = f(x + k2x / 2.0, v + k2v / 2.0, t + DELTA_T / 2.0) * DELTA_T / M;
  let k3x = (v + k2v / 2.0) * DELTA_T;
  let k4v = f(x + k3x, v + k3v, t + DELTA_T) * DELTA_T / M;
  let k4x = (v + k3v) * DELTA_T;
  let v = (k1v + 2.0 * k2v + 2.0 * k3v + k4v) / 6.0;
  let x = (k1x + 2.0 * k2x + 2.0 * k3x + k4x) / 6.0;
  (x, v)
}

fn f(x: f64, _v: f64, _t: f64) -> f64 {
  - K * x
}

fn plot(x: Vec<f64>, y: Vec<f64>) {
  let mut fg = Figure::new();
  fg.set_terminal("pngcairo", "output.png");
  fg.axes2d()
      .lines(&x, &y, &[Caption("A line"), Color("black")]);
  fg.show();
}

【Rust】Rustで数値計算プロジェクトを試す その5

昨日の件を友達に話したところ、crate側がRustのバージョンに対応してなかった可能性が高いようです。

これってちょっと見てみたけど、crate側の実装の問題な気がする。gnuplotのラッパーの方使ったら出力できた。

https://www.mathgram.xyz/entry/rust/gd

去年の中盤から更新されてなかったから

https://www.utam0k.jp/blog/2018/05/28/rust_std_default/

この一番下のエラーな気がするんだよな

というわけで、次のgnuplotライブラリを使って、サインカーブを出力してみます。

docs.rs

use gnuplot::{Figure, Caption, Color};

fn main() {
  let mut x = vec![];
  let mut y = vec![];

  for i in 0..1000 {
    let v = (i as f64) / 100.0;
    x.push(v);
    y.push(v.sin());
  }
  plot(x, y);
}

fn plot(x: Vec<f64>, y: Vec<f64>) {
  let mut fg = Figure::new();
  fg.set_terminal("pngcairo", "test.png");
  fg.axes2d()
      .lines(&x, &y, &[Caption("A line"), Color("black")]);
  fg.show();
}

できました！これで「サインカーブをプロットする」のというトイケースはクリアですね。

【Rust】Rustで数値計算プロジェクトを試す その4

昨日の続きで、Rustのcratesを利用してサインカーブをプロットしたいのですが、エラーが出てしまっています。gitconfigは、以前source treeを利用していた際のものが残っていただけだったので、ちゃんと設定し直しました。

エラー箇所を確認するためのシンプルな例をコンパイルしてみました。

extern crate plotlib;
use plotlib::scatter::Scatter;

fn main() {
    // Scatter plots expect a list of pairs
    let data1 = [(-3.0, 2.3), (-1.6, 5.3), (0.3, 0.7), (4.3, -1.4), (6.4, 4.3), (8.5, 3.7)];

    // We create our scatter plot from the data
    let s1 = Scatter::from_slice(&data1);
    println!("{:?}", s1);
}

やはりこういうエラーが出ます。

$ cargo run --verbose
       Fresh svg v0.5.12
       Fresh plotlib v0.3.0
   Compiling duplex-pendulum v0.1.0 (/Users/takeshi/Desktop/projects/duplex-pendulum)
     Running `rustc --edition=2018 --crate-name duplex_pendulum src/main.rs --color always --crate-type bin --emit=dep-info,link -C debuginfo=2 -C metadata=f4638e3323e9a474 -C extra-filename=-f4638e3323e9a474 --out-dir /Users/takeshi/Desktop/projects/duplex-pendulum/target/debug/deps -C incremental=/Users/takeshi/Desktop/projects/duplex-pendulum/target/debug/incremental -L dependency=/Users/takeshi/Desktop/projects/duplex-pendulum/target/debug/deps--extern plotlib=/Users/takeshi/Desktop/projects/duplex-pendulum/target/debug/deps/libplotlib-2ce36c14f3ccea1e.rlib`
error[E0599]: no function or associated item named `from_slice` found for type `plotlib::scatter::Scatter` in the current scope
 --> src/main.rs:9:14
  |
9 |     let s1 = Scatter::from_slice(&data1);
  |              ^^^^^^^^^^^^^^^^^^^ function or associated item not found in `plotlib::scatter::Scatter`

error: aborting due to previous error

For more information about this error, try `rustc --explain E0599`.
error: Could not compile `duplex-pendulum`.

Caused by:
  process didn't exit successfully: `rustc --edition=2018 --crate-name duplex_pendulum src/main.rs --color always --crate-type bin --emit=dep-info,link -C debuginfo=2 -C metadata=f4638e3323e9a474 -C extra-filename=-f4638e3323e9a474 --out-dir /Users/takeshi/Desktop/projects/duplex-pendulum/target/debug/deps -C incremental=/Users/takeshi/Desktop/projects/duplex-pendulum/target/debug/incremental -L dependency=/Users/takeshi/Desktop/projects/duplex-pendulum/target/debug/deps --extern plotlib=/Users/takeshi/Desktop/projects/duplex-pendulum/target/debug/deps/libplotlib-2ce36c14f3ccea1e.rlib` (exit code: 1)

ただ、他の関数（style）を無理やり呼び出すコードに変更すると「function or associated item not found in〜」のエラーは起きないし（引数の型の問題でエラー起きる）、名前空間関連の問題であるように思います。

github.com

なんだか消化不良ですが今日はこれくらいで。Rustの最近のコンパイラのバージョンで名前空間の仕様が変わったとかかな…。

【Rust】Rustで数値計算プロジェクトを試す その3

数値計算のプロジェクトでは、結局グラフの形式で出力することも多いと思います。csvで出力して古き良きgnuplotでもいいんですが、Pythonのmatplotlibみたいなのが無いか調べます。

pure Rustのこちらのライブラリを試してみます。

github.com

しかし、こちらのサンプルコードをコンパイルしようとするとエラーが出てしまいます。

extern crate plotlib;
use plotlib::scatter::Scatter;
use plotlib::scatter;
use plotlib::style::{Marker, Point};
use plotlib::view::View;
use plotlib::page::Page;

fn main() {
    // Scatter plots expect a list of pairs
    let data1 = [(-3.0, 2.3), (-1.6, 5.3), (0.3, 0.7), (4.3, -1.4), (6.4, 4.3), (8.5, 3.7)];

    // We create our scatter plot from the data
    let s1 = Scatter::from_slice(&data1)
        .style(scatter::Style::new()
            .marker(Marker::Square) // setting the marker to be a square
            .colour("#DD3355")); // and a custom colour

    // We can plot multiple data sets in the same view
    let data2 = [(-1.4, 2.5), (7.2, -0.3)];
    let s2 = Scatter::from_slice(&data2)
        .style(scatter::Style::new() // uses the default marker
            .colour("#35C788")); // and a different colour

    // The 'view' describes what set of data is drawn
    let v = View::new()
        .add(&s1)
        .add(&s2)
        .x_range(-5., 10.)
        .y_range(-2., 6.)
        .x_label("Some varying variable")
        .y_label("The response of something");

    // A page with a single view is then saved to an SVG file
    Page::single(&v).save("scatter.svg");
}

いくつかエラー出てますが、「どうやら関数が名前空間に存在しない」ような感じのエラーのようです。ただ、実際のコード見てもあるんだよなあ…。

> cargo run
   Compiling duplex-pendulum v0.1.0 (/Users/takeshi/Desktop/projects/duplex-pendulum)
error[E0599]: no function or associated item named `from_slice` found for type `plotlib::scatter::Scatter` in the current scope
  --> src/main.rs:13:14
   |
13 |     let s1 = Scatter::from_slice(&data1)
   |              ^^^^^^^^^^^^^^^^^^^ function or associated item not found in `plotlib::scatter::Scatter`

また、私も ~/.gitconfig の設定の問題で、↓の問題が起こってcratesを追加できなくなってしまっていました。これは解決したのですが、いろいろ設定いじってたらgithubの認証が利用できなくなる問題が…。後でなんとかしよう。

github.com

【Rust】Rustで数値計算プロジェクトを試す その2

前回のあらすじ: 「Rustで数値計算（二重振り子）して結果をアニメーションで表示したい！」と思っています。クソ記事量産してる感もあるんですが、以前チュートリアルをしていて「どこまで自分が分かっていて、どこから先分かってないのか分からない」みたいな状況になってしまったので振り返りのために書いてます。

簡単なトイケースとして「サインカーブをプロットする」のを目標にコードを書いています。ひとまずこんなの書けました。

fn main() {
    let mut curve = vec![];

    for i in 0..100 {
        curve.push((i as f64).sin());
    }

    println!("{:?}", curve);
}

ちょっと詰まったのが、 int のままでは sin() が呼び出せず、 float 型にキャストする必要があってその方法くらいですね。このあたりはC++の反省を活かして、型変換を明示する仕様に作り上げてるんだと勝手に思ってます。

あと、不要なライブラリをインポートしないイメージがあるんですが、 #include <math.h> みたいなこと書かなくてもいいってことは、余計なものをロードしない仕組みもありそうな雰囲気してます（雰囲気で話してます）。

「Rustにmatplotlibのようなライブラリが無いか」も探しています。

github.com

ただ、内部でPythonを呼び出してるライブラリが多く、「普段pipenv（pyenv + venv）の仮想環境使ってるけど、その呼び出しするときにけっこう詰まりそうだな」って思ってしまっています。

【Rust】Rustで数値計算プロジェクトを試す その1

これから伸びていって「今はPython(numpy)でやっているような数値計算もRustで、高速に抽象的に書けるようになるんじゃないかな」と思っています。

github.com

今まで何度かトライしたものの、「お前まだチュートリアルやってるのかよ」ってプログラミング友達に馬鹿にされる毎日でした。ひとまず使い方として想定する「数値計算プログラム」をRustで書いてみて、モチベーション保ちつつ実装していく作戦です。

Haskellで二重振り子のシミュレーションを書いた。実はこの振り子50本あって初期値を1兆分の1ぐらいずらしてあるので途中からめちゃくちゃ分岐する pic.twitter.com/3dh4Metcpj

— lotz△ (@lotz84_) 2017年5月22日

Haskellと同様にtraitで数学の抽象的な概念を扱えるらしいし、C++の代替言語を目指しているだけあってGC等のプログラマーがコントロールできない遅い処理が排除されています（その分lifetimeの管理が難しいのですが…）。

今日はひとまず、「配列を破壊的変更してprintする」ことをしました。初歩的ですが、ひとまず第一歩です。

$ rustup --version
rustup 1.16.0 (beab5ac2b 2018-12-06)

コードはこんな感じ。

fn main() {
    let mut v: Vec<i32> = vec![1, 2, 3];

    println!("v = {:?}", v);

    update_vector(&mut v);
    // debug print
    println!("v = {:?}", v);
}

fn update_vector(v: &mut Vec<i32>) {
    v.push(1);
    v[0] = 4;
    println!("v = {:?}", v);
}