Skip to content

Commit fdd85f9

Browse files
jackyarndleyjackyarndley
committed
further optimization of series approximation and packing
1 parent a53a4c1 commit fdd85f9

File tree

7 files changed

+117
-148
lines changed

7 files changed

+117
-148
lines changed

‎.gitignore‎

Lines changed: 4 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,7 @@
11
/target
22
**/*.rs.bk
33
.idea
4-
*.kfr
4+
*.kfr
5+
*.svg
6+
*.data
7+
*.old

‎Cargo.lock‎

Lines changed: 57 additions & 57 deletions
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

‎Cargo.toml‎

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,7 @@ rug = "1.9.0"
1313
itertools = "^0.9.0"
1414

1515
#Additional commands that can improve performance (maybe by around 5-10%)
16-
#[profile.release]
16+
[profile.release]
1717
#lto = "fat"
1818
#codegen-units = 1
19+
debug = true

‎output.png‎

-631 KB
Loading

‎src/math/perturbation.rs‎

Lines changed: 8 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
use crate::util::PixelData;
1+
use crate::util::{PixelData, FloatExp};
22

33
use rayon::prelude::*;
44
use crate::math::reference::Reference;
@@ -7,12 +7,12 @@ pub struct Perturbation {}
77

88
impl Perturbation {
99
pub fn iterate(pixel_data: &mut Vec<PixelData>, reference: &Reference, reference_current_iteration: usize) {
10-
pixel_data.par_chunks_mut(4)
10+
pixel_data.par_chunks_mut(64)
1111
.for_each(|pixel_data| {
1212
for pixel in pixel_data {
1313
let mut scaled_iterations = 0;
14-
let mut scaled_scale_factor_1 = 2.0f64.powi(pixel.delta_current.exponent);
15-
let mut scaled_delta_reference = 2.0f64.powi(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
14+
let mut scaled_scale_factor_1 = 1.0f64.ldexp(pixel.delta_current.exponent);
15+
let mut scaled_delta_reference = 1.0f64.ldexp(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
1616

1717
while pixel.iteration < reference_current_iteration {
1818
let delta_current_float = scaled_scale_factor_1 * pixel.delta_current.mantissa;
@@ -43,8 +43,8 @@ impl Perturbation {
4343
// reset the scaled counter
4444
pixel.delta_current.reduce();
4545

46-
scaled_scale_factor_1 = 2.0f64.powi(pixel.delta_current.exponent);
47-
scaled_delta_reference = 2.0f64.powi(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
46+
scaled_scale_factor_1 = 1.0f64.ldexp(pixel.delta_current.exponent);
47+
scaled_delta_reference = 1.0f64.ldexp(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
4848

4949
scaled_iterations = 0;
5050
},
@@ -61,8 +61,8 @@ impl Perturbation {
6161
if scaled_iterations > 250 {
6262
pixel.delta_current.reduce();
6363

64-
scaled_scale_factor_1 = 2.0f64.powi(pixel.delta_current.exponent);
65-
scaled_delta_reference = 2.0f64.powi(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
64+
scaled_scale_factor_1 = 1.0f64.ldexp(pixel.delta_current.exponent);
65+
scaled_delta_reference = 1.0f64.ldexp(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
6666

6767
scaled_iterations = 0;
6868
}

‎src/renderer.rs‎

Lines changed: 26 additions & 25 deletions
Original file line numberDiff line numberDiff line change
@@ -102,27 +102,29 @@ impl FractalRenderer {
102102

103103
let time = Instant::now();
104104

105-
let indices = (0..self.image_width).cartesian_product(0..self.image_height);
106-
107-
let mut pixel_data = indices.into_iter().par_bridge()
108-
.map(|(i, j)| {
109-
let element = ComplexFixed::new(i as f64 * delta_pixel + delta_top_left.re, j as f64 * delta_pixel + delta_top_left.im);
110-
let point_delta = ComplexExtended::new(element, -self.zoom.exponent);
111-
let new_delta = series_approximation.evaluate(point_delta);
112-
113-
PixelData {
114-
image_x: i,
115-
image_y: j,
116-
iteration: reference.start_iteration,
117-
delta_centre: point_delta,
118-
delta_reference: point_delta,
119-
delta_start: new_delta,
120-
delta_current: new_delta,
121-
derivative_current: ComplexFixed::new(1.0, 0.0),
122-
glitched: false,
123-
escaped: false
124-
}
125-
}).collect::<Vec<PixelData>>();
105+
let indices = (0..self.image_width).cartesian_product(0..self.image_height).collect::<Vec<(usize, usize)>>();
106+
107+
let mut pixel_data = Vec::with_capacity(self.image_width * self.image_height);
108+
109+
pixel_data = indices.into_par_iter()
110+
.map(|(i, j)| {
111+
let element = ComplexFixed::new(i as f64 * delta_pixel + delta_top_left.re, j as f64 * delta_pixel + delta_top_left.im);
112+
let point_delta = ComplexExtended::new(element, -self.zoom.exponent);
113+
let new_delta = series_approximation.evaluate(point_delta);
114+
115+
PixelData {
116+
image_x: i,
117+
image_y: j,
118+
iteration: reference.start_iteration,
119+
delta_centre: point_delta,
120+
delta_reference: point_delta,
121+
delta_start: new_delta,
122+
delta_current: new_delta,
123+
derivative_current: ComplexFixed::new(1.0, 0.0),
124+
glitched: false,
125+
escaped: false
126+
}
127+
}).collect::<Vec<PixelData>>();
126128

127129
println!("{:<14}{:>6} ms", "Packing", time.elapsed().as_millis());
128130

@@ -131,7 +133,7 @@ impl FractalRenderer {
131133
println!("{:<14}{:>6} ms", "Iteration", time.elapsed().as_millis());
132134

133135
let time = Instant::now();
134-
Colouring::IterationSmooth.run(&pixel_data, &mut self.image, self.maximum_iteration, delta_pixel, &reference);
136+
Colouring::Iteration.run(&pixel_data, &mut self.image, self.maximum_iteration, delta_pixel, &reference);
135137
println!("{:<14}{:>6} ms", "Coloring", time.elapsed().as_millis());
136138

137139
let time = Instant::now();
@@ -160,17 +162,16 @@ impl FractalRenderer {
160162
for data in pixel_data {
161163
data.iteration = reference.start_iteration;
162164
data.glitched = false;
163-
data.escaped = false;
164165
data.delta_current = data.delta_start - delta_z;
165166
data.delta_reference = data.delta_centre - reference_wrt_sa;
166167
// might not need the evaluate here as if we store it separately, there is no need
167-
data.derivative_current = ComplexFixed::new(1.0, 0.0);
168+
// data.derivative_current = ComplexFixed::new(1.0, 0.0);
168169
}
169170
});
170171

171172
Perturbation::iterate(&mut pixel_data, &r, r.current_iteration);
172173

173-
Colouring::IterationSmooth.run(&pixel_data, &mut self.image, self.maximum_iteration, delta_pixel, &r);
174+
Colouring::Iteration.run(&pixel_data, &mut self.image, self.maximum_iteration, delta_pixel, &r);
174175

175176
// Remove all non-glitched points from the remaining points
176177
pixel_data.retain(|packet| {

‎src/util/complex_extended.rs‎

Lines changed: 20 additions & 56 deletions
Original file line numberDiff line numberDiff line change
@@ -14,24 +14,18 @@ pub struct ComplexExtended {
1414
impl ComplexExtended {
1515
#[inline]
1616
pub fn new(mantissa: Complex<f64>, exponent: i32) -> Self {
17-
let mut temp = ComplexExtended {
17+
ComplexExtended {
1818
mantissa,
1919
exponent
20-
};
21-
22-
temp.reduce();
23-
temp
20+
}
2421
}
2522

2623
#[inline]
2724
pub fn new2(re: f64, im: f64, exponent: i32) -> Self {
28-
let mut temp = ComplexExtended {
25+
ComplexExtended {
2926
mantissa: Complex::<f64>::new(re, im),
3027
exponent
31-
};
32-
33-
temp.reduce();
34-
temp
28+
}
3529
}
3630

3731
#[inline]
@@ -49,7 +43,7 @@ impl ComplexExtended {
4943

5044
#[inline]
5145
pub fn to_float(&self) -> Complex<f64> {
52-
self.mantissa * 2.0f64.powi(self.exponent)
46+
self.mantissa * 1.0f64.ldexp(self.exponent)
5347
}
5448

5549
#[inline]
@@ -71,36 +65,25 @@ impl ComplexExtended {
7165
impl AddAssign for ComplexExtended {
7266
#[inline]
7367
fn add_assign(&mut self, other: Self) {
74-
if self.mantissa.re == 0.0 && self.mantissa.im == 0.0 {
75-
self.mantissa = other.mantissa;
76-
self.exponent = other.exponent;
77-
} else if other.mantissa.re == 0.0 && other.mantissa.im == 0.0 {
78-
return;
68+
if self.exponent > other.exponent {
69+
self.mantissa += other.mantissa * 1.0f64.ldexp(other.exponent - self.exponent);
7970
} else {
80-
if self.exponent == other.exponent {
81-
self.mantissa += other.mantissa;
82-
} else if self.exponent > other.exponent {
83-
self.mantissa += other.mantissa / 2.0f64.powi(self.exponent - other.exponent);
84-
} else {
85-
self.mantissa /= 2.0f64.powi(other.exponent - self.exponent);
86-
self.exponent = other.exponent;
87-
self.mantissa += other.mantissa;
88-
}
71+
self.mantissa *= 1.0f64.ldexp(self.exponent - other.exponent);
72+
self.mantissa += other.mantissa;
73+
self.exponent = other.exponent;
8974
}
9075
}
9176
}
9277

9378
impl SubAssign for ComplexExtended {
9479
#[inline]
9580
fn sub_assign(&mut self, other: Self) {
96-
if self.exponent == other.exponent {
97-
self.mantissa -= other.mantissa;
98-
} else if self.exponent > other.exponent {
99-
self.mantissa -= other.mantissa / 2.0f64.powi(self.exponent - other.exponent);
81+
if self.exponent > other.exponent {
82+
self.mantissa -= other.mantissa * 1.0f64.ldexp(other.exponent - self.exponent);
10083
} else {
101-
self.mantissa /= 2.0f64.powi(other.exponent - self.exponent);
102-
self.exponent = other.exponent;
84+
self.mantissa *= 1.0f64.ldexp(self.exponent - other.exponent);
10385
self.mantissa -= other.mantissa;
86+
self.exponent = other.exponent;
10487
}
10588
}
10689
}
@@ -125,7 +108,6 @@ impl DivAssign for ComplexExtended {
125108
fn div_assign(&mut self, other: Self) {
126109
self.mantissa /= other.mantissa;
127110
self.exponent -= other.exponent;
128-
self.reduce();
129111
}
130112
}
131113

@@ -134,19 +116,10 @@ impl Add<ComplexExtended> for ComplexExtended {
134116

135117
#[inline]
136118
fn add(self, other: Self) -> Self::Output {
137-
if self.mantissa.re == 0.0 && self.mantissa.im == 0.0 {
138-
other
139-
} else if other.mantissa.re == 0.0 && other.mantissa.im == 0.0 {
140-
self
119+
if self.exponent > other.exponent {
120+
ComplexExtended::new(self.mantissa + other.mantissa * 1.0f64.ldexp(other.exponent - self.exponent), self.exponent)
141121
} else {
142-
let (new_mantissa, new_exponent) = if self.exponent == other.exponent {
143-
(self.mantissa + other.mantissa, self.exponent)
144-
} else if self.exponent > other.exponent {
145-
(self.mantissa + other.mantissa / 2.0f64.powi(self.exponent - other.exponent), self.exponent)
146-
} else {
147-
(other.mantissa + self.mantissa / 2.0f64.powi(other.exponent - self.exponent), other.exponent)
148-
};
149-
ComplexExtended::new(new_mantissa, new_exponent)
122+
ComplexExtended::new(other.mantissa + self.mantissa * 1.0f64.ldexp(self.exponent - other.exponent), other.exponent)
150123
}
151124
}
152125
}
@@ -156,19 +129,10 @@ impl Sub<ComplexExtended> for ComplexExtended {
156129

157130
#[inline]
158131
fn sub(self, other: Self) -> Self::Output {
159-
if self.mantissa.re == 0.0 && self.mantissa.im == 0.0 {
160-
ComplexExtended::new(-other.mantissa, other.exponent)
161-
} else if other.mantissa.re == 0.0 && other.mantissa.im == 0.0 {
162-
self
132+
if self.exponent > other.exponent {
133+
ComplexExtended::new(self.mantissa - other.mantissa * 1.0f64.ldexp(other.exponent - self.exponent), self.exponent)
163134
} else {
164-
let (new_mantissa, new_exponent) = if self.exponent == other.exponent {
165-
(self.mantissa - other.mantissa, self.exponent)
166-
} else if self.exponent > other.exponent {
167-
(self.mantissa - other.mantissa / 2.0f64.powi(self.exponent - other.exponent), self.exponent)
168-
} else {
169-
(-1.0 * other.mantissa + self.mantissa / 2.0f64.powi(other.exponent - self.exponent), other.exponent)
170-
};
171-
ComplexExtended::new(new_mantissa, new_exponent)
135+
ComplexExtended::new(other.mantissa - self.mantissa * 1.0f64.ldexp(self.exponent - other.exponent), other.exponent)
172136
}
173137
}
174138
}

0 commit comments

Comments
 (0)