LoginSignup
6

More than 3 years have passed since last update.

Rustで数値計算してみた

Last updated at Posted at 2019-06-05

RustはSIMDをサポートしているため数値計算も視野に入れていると見ました。現在SIMDを試せる環境がないのですが,システム制御用の信号処理器(とりあえず微分と積分を実現するIIRフィルタ)を作ってみました。

main.rs

use std::io::prelude::*;
use std::fs::File;
use std::io::BufWriter;
use std::f32;

mod integrator1;
mod differentiator1;

use integrator1::Integrator1;
use differentiator1::Differentiator1;

fn main() { 

    let ts=1e-4;
    let mut t:f32=0.0;

    let sim_range=20000;

    let fp1 = format!("dif.dat");
    let fp2 = format!("int.dat");

    let mut buffer1 = File::create(fp1).unwrap();
    let mut buffer2 = File::create(fp2).unwrap();

    //Buffering
    let mut buffer1 = BufWriter::new(buffer1);
    let mut buffer2 = BufWriter::new(buffer2);

    let mut int=Integrator1::new(ts);
    let mut dif=Differentiator1::new(ts,100.0);

    let freq:f32=5.0*f32::consts::PI;
    let mut uix:f32;    //input for integrator
    let mut udx:f32;    //input for differentiator
    let mut ix:f32;     //output of integrator
    let mut dx:f32;     //output of differentiator

    for _ in 0..sim_range {
        uix=(freq*t).sin();
        udx=1.0-(freq*t).cos();

        ix=int.get(uix);
        dx=dif.get(udx);

        write!(buffer1, "{:.04e} {:.04e} {:.04e}\n",t,udx,dx).unwrap();
        write!(buffer2, "{:.04e} {:.04e} {:.04e}\n",t,uix,ix).unwrap();

        t+=ts;
    }
}

IIRフィルタ実装にはstructを利用してみました。これで複数のIIRを呼べるようになりました。

differentiaotr1.rs
use std::fmt;

pub struct Differentiator1{
    yz1:f32,
    uz1:f32,
    a:f32,
    b:f32,
}

impl fmt::Display for Differentiator1 {
    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
        write!(f, "{:.06e}", self.yz1)
    }
}

impl Differentiator1{
    pub fn new(ts:f32,g:f32)->Differentiator1{
        Differentiator1{
            yz1:0.0,
            uz1:0.0,
            a:2.0*g/(2.0+g*ts),
            b:(2.0-g*ts)/(2.0+g*ts),
        }
    }

    pub fn _zero(&self)->Differentiator1{
        Differentiator1{
            yz1:0.0,
            uz1:0.0,
            a:self.a,
            b:self.b,
        }
    }

    pub fn get(&mut self,u:f32)->f32{
        self.yz1=(u-self.uz1)*self.a+self.yz1*self.b;
        self.uz1=u;
        self.yz1
    }
}
integrator1.rs
use std::fmt;

pub struct Integrator1{
    yz1:f32,
    uz1:f32,
    a:f32,
}

impl fmt::Display for Integrator1 {
    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
        write!(f, "{:.06e}", self.yz1)
    }
}

impl Integrator1{
    pub fn new(ts:f32)->Integrator1{
        Integrator1{
            yz1:0.0,
            uz1:0.0,
            a:ts/2.0,
        }
    }

    pub fn _zero(&self)->Integrator1{
        Integrator1{
            yz1:0.0,
            uz1:0.0,
            a:self.a
        }
    }

    pub fn get(&mut self,u:f32)->f32{
        self.yz1=(u+self.uz1)*self.a+self.yz1;
        self.uz1=u;
        self.yz1
    }
}

入出力結果を見ると,微積分ができていました。
xx_res.png
c++より簡単に実装できる気がしました。

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
6