ising_model/
ising_model.rs

1use rapl::utils::random::NdarrRand;
2use rapl::*;
3use std::{
4    io::{stdout, Write},
5    thread::sleep,
6    time::Duration,
7};
8
9const N: usize = 18; // Number of spins in each direction.
10const T: f32 = 0.05; // Temperature
11const STEPS: usize = 700; //Steps of simulation
12const M: usize = 20; //Number of updates per metropolis circle
13
14fn metropolis(spin_arr: &mut Ndarr<f32, U2>) {
15    let energy: Ndarr<f32, U2> = 2.
16        * spin_arr.clone()
17        * (spin_arr.roll(1, 0) + spin_arr.roll(-1, 0) + spin_arr.roll(1, 1) + spin_arr.roll(-1, 1));
18    let temp_exp = (-&energy / T).exp();
19    let indexes: Vec<usize> = (0..N).collect();
20    let i_s = NdarrRand::choose(&indexes, [M], None); //random i indexes
21    let j_s = NdarrRand::choose(&indexes, [M], None); //random j indexes
22    let p_swich: f32 = NdarrRand::uniform(0.0, 1.0, [1], None)[0]; // selection probability of random spin switch
23    for (i, j) in i_s.data.iter().zip(j_s.data.iter()) {
24        if energy[[*i, *j]] < 0.0 || p_swich < temp_exp[[*i, *j]] {
25            spin_arr[[*i, *j]] *= &-1.0;
26        }
27    }
28}
29
30fn main() {
31    let mut stdout = stdout();
32    let mut spin_arr = NdarrRand::choose(&[-1.0, 1.0], [N, N], None);
33    stdout.flush().unwrap();
34    stdout.write_all(b"\x1B[2J\x1B[1;1H").unwrap();
35    for i in 0..STEPS {
36        metropolis(&mut spin_arr); //updates array with metropolis algorithm
37        let vis = spin_arr.map(|x| {
38            if *x < 0.0 {
39                "░".to_string()
40            } else {
41                "█".to_string()
42            }
43        }); //make it pretty
44
45        println!("{}", vis);
46
47        println!(
48            "\n The Ising Model using rapl: [Step {} out of {}] \n \n for more info: https://en.wikipedia.org/wiki/Ising_model"
49        , i + 1, STEPS);
50        sleep(Duration::from_millis(5));
51        stdout.write_all(b"\x1B[1;1H").unwrap();
52        stdout.flush().unwrap();
53    }
54    stdout.write_all(b"\x1B[2J\x1B[1;1H").unwrap();
55    stdout.flush().unwrap();
56}