Skip to content

Commit d34a350

Browse files
jackyarndleyjackyarndley
committed
implemented analytic derivatives
1 parent 88098ae commit d34a350

File tree

7 files changed

+194
-90
lines changed

7 files changed

+194
-90
lines changed

‎high.toml‎

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -20,4 +20,6 @@ data_storage_interval = 100
2020
valid_iteration_frame_multiplier = 0.10
2121
valid_iteration_probe_multiplier = 0.01
2222

23-
experimental = true
23+
experimental = true
24+
25+
analytic_derivative = true

‎locations/e10000.toml‎

Lines changed: 1 addition & 1 deletion
Large diffs are not rendered by default.

‎src/math/perturbation.rs‎

Lines changed: 73 additions & 48 deletions
Original file line numberDiff line numberDiff line change
@@ -6,7 +6,7 @@ use crate::math::reference::Reference;
66
pub struct Perturbation {}
77

88
impl Perturbation {
9-
pub fn iterate(pixel_data: &mut [PixelData], reference: &Reference) {
9+
pub fn iterate_normal(pixel_data: &mut [PixelData], reference: &Reference) {
1010
pixel_data.par_chunks_mut(4)
1111
.for_each(|pixel_data| {
1212
for pixel in pixel_data {
@@ -58,8 +58,6 @@ impl Perturbation {
5858
} else {
5959
// If the reference is not small, use the usual method
6060

61-
// pixel.delta_current.mantissa = 2.0 * reference.data[pixel.iteration - reference.start_iteration].z_fixed * pixel.delta_current.mantissa + temp * pixel.delta_current.mantissa + scaled_delta_reference;
62-
6361
// 4 multiplications and 2 additions
6462
pixel.delta_current.mantissa *= z + reference_data.z_fixed;
6563
// 2 additions
@@ -82,67 +80,94 @@ impl Perturbation {
8280
if !pixel.escaped && !pixel.glitched {
8381
pixel.iteration = reference.current_iteration;
8482
}
83+
}
84+
});
85+
}
86+
87+
pub fn iterate_normal_plus_derivative(pixel_data: &mut [PixelData], reference: &Reference) {
88+
pixel_data.par_chunks_mut(4)
89+
.for_each(|pixel_data| {
90+
for pixel in pixel_data {
91+
let mut scaled_iterations = 0;
92+
let mut scaled_scale_factor_1 = 1.0f64.ldexp(pixel.delta_current.exponent);
93+
let mut scaled_scale_factor_2 = 1.0f64.ldexp(-pixel.derivative_current.exponent);
94+
let mut scaled_delta_reference = 1.0f64.ldexp(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
95+
96+
let val1 = pixel.iteration - reference.start_iteration;
97+
let val2 = reference.current_iteration - reference.start_iteration;
98+
let reference_slice = &reference.reference_data[val1..=val2];
8599

86-
// while pixel.iteration < reference.current_iteration {
87-
// let reference_data = &reference.reference_data[pixel.iteration - reference.start_iteration];
100+
for additional_iterations in 0..=(reference.current_iteration - pixel.iteration) {
101+
let reference_data = &reference_slice[additional_iterations];
88102

89-
// // 2 multiplications and 2 adds
90-
// let z = reference_data.z_fixed + scaled_scale_factor_1 * pixel.delta_current.mantissa;
103+
// 2 multiplications and 2 adds
104+
let z = reference_data.z_fixed + scaled_scale_factor_1 * pixel.delta_current.mantissa;
91105

92-
// if pixel.delta_current.exponent > -500 {
93-
// // 2 multiplications and one add
94-
// let z_norm = z.norm_sqr();
106+
if pixel.delta_current.exponent > -500 {
107+
// 2 multiplications and one add
108+
let z_norm = z.norm_sqr();
95109

96-
// if z_norm < reference_data.z_tolerance {
97-
// pixel.glitched = true;
98-
// break;
99-
// }
110+
if z_norm < reference_data.z_tolerance {
111+
pixel.iteration += additional_iterations;
112+
pixel.glitched = true;
113+
break;
114+
}
100115

101-
// if z_norm > 1e16 {
102-
// pixel.escaped = true;
103-
// pixel.delta_current.mantissa = pixel.delta_current.to_float();
104-
// pixel.delta_current.exponent = 0;
105-
// break;
106-
// }
107-
// }
116+
if z_norm > 1e16 {
117+
pixel.iteration += additional_iterations;
118+
pixel.escaped = true;
119+
pixel.delta_current.mantissa = pixel.delta_current.to_float();
120+
pixel.delta_current.exponent = 0;
121+
break;
122+
}
123+
124+
pixel.derivative_current.mantissa *= 2.0 * z;
125+
pixel.derivative_current.mantissa += scaled_scale_factor_2;
126+
}
108127

109-
// if reference_data.extended_precision_required {
110-
// // If the reference is small, use the slow extended method
111-
// pixel.delta_current *= reference_data.z_extended * 2.0 + pixel.delta_current;
112-
// pixel.delta_current += pixel.delta_reference;
128+
if reference_data.extended_precision_required {
129+
// If the reference is small, use the slow extended method
130+
pixel.delta_current *= reference_data.z_extended * 2.0 + pixel.delta_current;
131+
pixel.delta_current += pixel.delta_reference;
113132

114-
// // reset the scaled counter
115-
// pixel.delta_current.reduce();
133+
// reset the scaled counter
134+
pixel.delta_current.reduce();
135+
pixel.derivative_current.reduce();
116136

117-
// scaled_scale_factor_1 = 1.0f64.ldexp(pixel.delta_current.exponent);
118-
// scaled_delta_reference = 1.0f64.ldexp(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
137+
scaled_scale_factor_1 = 1.0f64.ldexp(pixel.delta_current.exponent);
138+
scaled_scale_factor_2 = 1.0f64.ldexp(-pixel.derivative_current.exponent);
139+
scaled_delta_reference = 1.0f64.ldexp(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
119140

120-
// scaled_iterations = 0;
121-
// } else {
122-
// // If the reference is not small, use the usual method
141+
scaled_iterations = 0;
142+
} else {
143+
// If the reference is not small, use the usual method
123144

124-
// // pixel.delta_current.mantissa = 2.0 * reference.data[pixel.iteration - reference.start_iteration].z_fixed * pixel.delta_current.mantissa + temp * pixel.delta_current.mantissa + scaled_delta_reference;
145+
// 4 multiplications and 2 additions
146+
pixel.delta_current.mantissa *= z + reference_data.z_fixed;
147+
// 2 additions
148+
pixel.delta_current.mantissa += scaled_delta_reference;
125149

126-
// // 4 multiplications and 2 additions
127-
// pixel.delta_current.mantissa *= z + reference_data.z_fixed;
128-
// // 2 additions
129-
// pixel.delta_current.mantissa += scaled_delta_reference;
150+
scaled_iterations += 1;
130151

131-
// scaled_iterations += 1;
152+
// check the counter, if it is > 250, do a normalisation
153+
if scaled_iterations > 250 {
154+
pixel.delta_current.reduce();
155+
pixel.derivative_current.reduce();
132156

133-
// // check the counter, if it is > 250, do a normalisation
134-
// if scaled_iterations > 250 {
135-
// pixel.delta_current.reduce();
157+
scaled_scale_factor_1 = 1.0f64.ldexp(pixel.delta_current.exponent);
158+
scaled_scale_factor_2 = 1.0f64.ldexp(-pixel.derivative_current.exponent);
159+
scaled_delta_reference = 1.0f64.ldexp(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
136160

137-
// scaled_scale_factor_1 = 1.0f64.ldexp(pixel.delta_current.exponent);
138-
// scaled_delta_reference = 1.0f64.ldexp(pixel.delta_reference.exponent - pixel.delta_current.exponent) * pixel.delta_reference.mantissa;
161+
scaled_iterations = 0;
162+
}
163+
}
164+
}
139165

140-
// scaled_iterations = 0;
141-
// }
142-
// }
166+
pixel.derivative_current.reduce();
143167

144-
// pixel.iteration += 1;
145-
// }
168+
if !pixel.escaped && !pixel.glitched {
169+
pixel.iteration = reference.current_iteration;
170+
}
146171
}
147172
});
148173
}

‎src/math/series_approximation.rs‎

Lines changed: 16 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -417,29 +417,24 @@ impl SeriesApproximation {
417417
approximation
418418
}
419419

420-
// pub fn evaluate_derivative(&self, point_delta: ComplexExtended) -> FloatExtended {
421-
// let mut original_point_derivative_n = ComplexExtended::new(1.0, 0, 0.0, 0);
422-
// let mut approximation_derivative = ComplexExtended::new(0.0, 0, 0.0, 0);
423-
424-
// for k in 1..=self.order {
425-
// approximation_derivative += k as f64 * self.coefficients[k] * original_point_derivative_n;
426-
// original_point_derivative_n *= ComplexExtended::new(point_delta.re, 0, point_delta.im, 0);
427-
// };
428-
429-
// approximation_derivative.to_float()
430-
431-
432-
// let mut approximation_derivative = self.coefficients[self.order] * self.order as f64;
420+
pub fn evaluate_derivative(&self, point_delta: ComplexExtended, iteration: usize) -> ComplexExtended {
421+
// This could be improved to use the iteration option better
422+
// this assumes that the requested iteration is a multiple of the data interval
433423

434-
// for k in (1..=(self.order - 1)).rev() {
435-
// approximation *= point_delta;
436-
// approximation += self.coefficients[k] * k as f64;
437-
// }
424+
// 101 -> 100 / 100 = 1, 1 -> 0 / 100 = 0, 201 -> 200 / 100 = 2
425+
let new_coefficients = &self.coefficients[(iteration - 1) / self.data_storage_interval];
426+
// Horner's rule
427+
let mut approximation = new_coefficients[self.order];
428+
approximation *= self.order as f64;
438429

439-
// approximation *= point_delta;
440-
// approximation.reduce();
441-
// approximation
430+
for k in (1..=(self.order - 1)).rev() {
431+
approximation *= point_delta;
432+
approximation += new_coefficients[k] * k as f64;
433+
}
442434

435+
// println!("{:?}", approximation);
443436

444-
// }
437+
// approximation *= point_delta;
438+
approximation
439+
}
445440
}

‎src/renderer.rs‎

Lines changed: 43 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -26,6 +26,7 @@ pub struct FractalRenderer {
2626
series_approximation: SeriesApproximation,
2727
render_indices: Vec<usize>,
2828
remove_centre: bool,
29+
analytic_derivative: bool,
2930
experimental: bool,
3031
}
3132

@@ -53,6 +54,8 @@ impl FractalRenderer {
5354
let valid_iteration_probe_multiplier = settings.get_float("valid_iteration_probe_multiplier").unwrap_or(0.02) as f32;
5455
let glitch_tolerance = settings.get_float("glitch_tolerance").unwrap_or(1.4e-6) as f64;
5556
let data_storage_interval = settings.get_int("data_storage_interval").unwrap_or(10) as usize;
57+
let analytic_derivative = settings.get_bool("analytic_derivative").unwrap_or(false);
58+
5659

5760
let data_type = match settings.get_str("export").unwrap_or(String::from("COLOUR")).to_ascii_uppercase().as_ref() {
5861
"RAW" | "EXR" => DataType::RAW,
@@ -122,7 +125,7 @@ impl FractalRenderer {
122125
auto_adjust_iterations,
123126
maximum_iteration,
124127
glitch_percentage,
125-
data_export: DataExport::new(image_width, image_height, display_glitches, data_type, palette, iteration_division),
128+
data_export: DataExport::new(image_width, image_height, display_glitches, data_type, palette, iteration_division, analytic_derivative),
126129
start_render_time: Instant::now(),
127130
remaining_frames,
128131
frame_offset,
@@ -131,6 +134,7 @@ impl FractalRenderer {
131134
series_approximation,
132135
render_indices,
133136
remove_centre,
137+
analytic_derivative,
134138
experimental
135139
}
136140
}
@@ -228,14 +232,20 @@ impl FractalRenderer {
228232
let point_delta = ComplexExtended::new(element, -self.zoom.exponent);
229233
let new_delta = self.series_approximation.evaluate(point_delta, chosen_iteration);
230234

235+
let derivative = if self.analytic_derivative {
236+
self.series_approximation.evaluate_derivative(point_delta, chosen_iteration)
237+
} else {
238+
ComplexExtended::new2(1.0, 0.0, 0)
239+
};
240+
231241
PixelData {
232242
image_x: i,
233243
image_y: j,
234244
iteration: chosen_iteration,
235245
delta_centre: point_delta,
236246
delta_reference: point_delta,
237247
delta_current: new_delta.clone(),
238-
derivative_current: ComplexFixed::new(1.0, 0.0),
248+
derivative_current: derivative,
239249
glitched: false,
240250
escaped: false,
241251
}
@@ -247,9 +257,13 @@ impl FractalRenderer {
247257
let iteration_time = Instant::now();
248258

249259
// This one has no offset because it is not a glitch resolving reference
250-
Perturbation::iterate(&mut pixel_data, &self.center_reference);
260+
if self.analytic_derivative {
261+
Perturbation::iterate_normal_plus_derivative(&mut pixel_data, &self.center_reference);
262+
} else {
263+
Perturbation::iterate_normal(&mut pixel_data, &self.center_reference);
264+
}
251265

252-
self.data_export.export_pixels(&pixel_data, self.maximum_iteration, &self.center_reference);
266+
self.data_export.export_pixels(&pixel_data, self.maximum_iteration, &self.center_reference, delta_pixel_extended);
253267
print!("| {:<15}", iteration_time.elapsed().as_millis());
254268
std::io::stdout().flush().unwrap();
255269

@@ -294,17 +308,33 @@ impl FractalRenderer {
294308
} else {
295309
let delta_current_reference = self.series_approximation.evaluate(glitch_reference_pixel.delta_centre, glitch_reference.start_iteration);
296310

297-
pixel_data.par_iter_mut()
298-
.for_each(|pixel| {
299-
pixel.iteration = glitch_reference.start_iteration;
300-
pixel.glitched = false;
301-
pixel.delta_current = self.series_approximation.evaluate( pixel.delta_centre, glitch_reference.start_iteration) - delta_current_reference;
302-
pixel.delta_reference = pixel.delta_centre - glitch_reference_pixel.delta_centre;
303-
});
311+
if self.analytic_derivative {
312+
pixel_data.par_iter_mut()
313+
.for_each(|pixel| {
314+
pixel.iteration = glitch_reference.start_iteration;
315+
pixel.glitched = false;
316+
pixel.delta_current = self.series_approximation.evaluate( pixel.delta_centre, glitch_reference.start_iteration) - delta_current_reference;
317+
pixel.delta_reference = pixel.delta_centre - glitch_reference_pixel.delta_centre;
318+
pixel.derivative_current = self.series_approximation.evaluate_derivative(pixel.delta_centre, glitch_reference.start_iteration);
319+
});
320+
} else {
321+
pixel_data.par_iter_mut()
322+
.for_each(|pixel| {
323+
pixel.iteration = glitch_reference.start_iteration;
324+
pixel.glitched = false;
325+
pixel.delta_current = self.series_approximation.evaluate( pixel.delta_centre, glitch_reference.start_iteration) - delta_current_reference;
326+
pixel.delta_reference = pixel.delta_centre - glitch_reference_pixel.delta_centre;
327+
});
328+
}
304329
};
305330

306-
Perturbation::iterate(&mut pixel_data, &glitch_reference);
307-
self.data_export.export_pixels(&pixel_data, self.maximum_iteration, &glitch_reference);
331+
if self.analytic_derivative {
332+
Perturbation::iterate_normal_plus_derivative(&mut pixel_data, &glitch_reference);
333+
} else {
334+
Perturbation::iterate_normal(&mut pixel_data, &glitch_reference);
335+
};
336+
337+
self.data_export.export_pixels(&pixel_data, self.maximum_iteration, &glitch_reference, delta_pixel_extended);
308338

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

0 commit comments

Comments
 (0)