1use std::io::prelude::*;
5
6extern crate ssimd;
7use ssimd::{u32x4, f32x4};
8
9#[inline(never)]
10fn mandelbrot_naive(c_x: f32, c_y: f32, max_iter: u32) -> u32 {
11 let mut x = c_x;
12 let mut y = c_y;
13 let mut count = 0;
14 while count < max_iter {
15 let xy = x * y;
16 let xx = x * x;
17 let yy = y * y;
18 let sum = xx + yy;
19 if sum > 4.0 {
20 break
21 }
22 count += 1;
23 x = xx - yy + c_x;
24 y = xy * 2.0 + c_y;
25 }
26 count
27}
28
29#[inline(never)]
30fn mandelbrot_vector(c_x: f32x4, c_y: f32x4, max_iter: u32) -> u32x4 {
31 let mut x = c_x;
32 let mut y = c_y;
33
34 let mut count = u32x4::splat(0);
35 for _ in 0..max_iter as usize {
36 let xy = x * y;
37 let xx = x * x;
38 let yy = y * y;
39 let sum = xx + yy;
40 let mask = sum.lt(f32x4::splat(4.0));
41
42 if !mask.any() { break }
43 count = count + mask.to_u().select(u32x4::splat(1),
44 u32x4::splat(0));
45
46 x = xx - yy + c_x;
47 y = xy + xy + c_y;
48 }
49 count
50}
51
52const COLOURS: &'static [(f32, f32, f32)] = &[(0.0, 7.0, 100.0),
53 (32.0, 107.0, 203.0),
54 (237.0, 255.0, 255.0),
55 (255.0, 170.0, 0.0),
56 (0.0, 2.0, 0.0)];
57const SCALE: f32 = 12.0;
58const LIMIT: u32 = 100;
59
60#[inline(never)]
61fn output_one(buf: &mut [u8], val: u32) {
62 let (r, g, b);
63 if val == LIMIT {
64 r = 0;
65 g = 0;
66 b = 0;
67 } else {
68 let val = (val as f32 % SCALE) * (COLOURS.len() as f32) / SCALE;
69 let left = val as usize % COLOURS.len();
70 let right = (left + 1) % COLOURS.len();
71
72 let p = val - left as f32;
73 let (r1, g1, b1) = COLOURS[left];
74 let (r2, g2, b2) = COLOURS[right];
75 r = (r1 + (r2 - r1) * p) as u8;
76 g = (g1 + (g2 - g1) * p) as u8;
77 b = (b1 + (b2 - b1) * p) as u8;
78 }
79 buf[0] = r;
80 buf[1] = g;
81 buf[2] = b;
82}
83
84fn main() {
85 let mut args = std::env::args();
86 args.next();
87 let width = args.next().unwrap().parse().unwrap();
88 let height = args.next().unwrap().parse().unwrap();
89
90 let left = -2.2;
91 let right = left + 3.0;
92 let top = 1.0;
93 let bottom = top - 2.0;
94
95 let width_step: f32 = (right - left) / width as f32;
96 let height_step: f32 = (bottom - top) / height as f32;
97
98 let adjust = f32x4::splat(width_step) * f32x4::new(0., 1., 2., 3.);
99
100 println!("P6 {} {} 255", width, height);
101 let mut line = vec![0; width * 3];
102
103 if args.next().is_none() {
104 for i in 0..height {
105 let y = f32x4::splat(top + height_step * i as f32);
106 for _j in 0..width/4 {
107 let j = 4*_j;
108 let x = f32x4::splat(left + width_step * j as f32) + adjust;
109 let ret = mandelbrot_vector(x, y, LIMIT);
110 for k in 0..4 { let val = ret.extract(k as u32); output_one(&mut line[3*(j + k)..3*(j + k + 1)], val); }
112 }
113 ::std::io::stdout().write(&line).unwrap();
114 }
115 } else {
116 for i in 0..height {
117 let y = top + height_step * i as f32;
118 for j in 0..width {
119 let x = left + width_step * j as f32;
120 let val = mandelbrot_naive(x, y, LIMIT);
121 output_one(&mut line[3*j..3*(j + 1)], val);
123 }
124 ::std::io::stdout().write(&line).unwrap();
125 }
126 }
127}