main.rs

#![feature(core_intrinsics)] //nightly needed

use std::intrinsics::{fadd_fast, fmul_fast};
use std::env;
use std::time::Instant;

fn main() {
    let args: Vec<String> = env::args().collect();

    let n: usize = args[1].trim().parse()
        .expect("Please type a number!");

    //arrays are immutable so I'll use Vec
    let mut a = vec![0.0 as f64; n * n];
    let mut b = vec![0.0 as f64; n * n];
    let mut c = vec![0.0 as f64; n * n];

    matrix_fill(n, &mut a);
    matrix_fill(n, &mut b);

    let t1 = Instant::now();
    matrix_mult(n, &a, &b, &mut c);
    let t2 = Instant::now();

    println!("{:?}\n{}", t2.duration_since(t1), c[n*(n/2) + n/2]);
}

fn matrix_mult(n: usize, a: &Vec<f64>, b: &Vec<f64>, c: &mut Vec<f64>) {
    let mut t = vec![0.0 as f64; n * n];

    for (i, b_row) in b.chunks_exact(n).enumerate() {
        for (j, &b_val) in b_row.iter().enumerate() {
            unsafe {
                *t.get_unchecked_mut(j*n + i) = b_val;
            }
        }
    }

    for (i, c_row) in c.chunks_exact_mut(n).enumerate() {
        for (j, c_val) in c_row.iter_mut().enumerate() {
            let a_row = &a[i*n..i*n+n];
            let t_row = &t[j*n..j*n+n];
            *c_val = a_row.iter().zip(t_row.iter()).fold(0.0, |sum, (&a,&t)| unsafe{ fadd_fast(fmul_fast(a,t), sum) } );
        }
    }

}

fn matrix_fill(n: usize, a: &mut Vec<f64>) {
    let tmp: f64 = 1.0 as f64 / n as f64 / n as f64;
    for i in 0..n {
        for j in 0..n {
	    a[i*n + j] = tmp * (i as i32 - j as i32) as f64 * (i + j) as f64;
	}
    }
}

fn matrix_print(n: usize, a: &Vec<f64>) {
    for i in 0..n {
        for j in 0..n {
	    print!("{:+07} ", a[i*n + j]);
	}
	println!();
    }
}