diff --git a/src/function.rs b/src/function.rs index 989aa0c..5226984 100644 --- a/src/function.rs +++ b/src/function.rs @@ -2,7 +2,7 @@ use crate::function_output::FunctionOutput; #[allow(unused_imports)] -use crate::misc::{debug_log, SteppedVector}; +use crate::misc::{debug_log, newtons_method, SteppedVector}; use crate::egui_app::{DEFAULT_FUNCION, DEFAULT_RIEMANN}; use crate::parsing::BackingFunction; @@ -228,109 +228,35 @@ impl FunctionEntry { // Finds roots fn roots(&mut self) { let resolution: f64 = (self.pixel_width as f64 / (self.max_x - self.min_x).abs()) as f64; - let mut root_list: Vec = Vec::new(); - let mut last_ele: Option = None; - for ele in self.output.back.as_ref().unwrap().iter() { - if last_ele.is_none() { - last_ele = Some(*ele); - continue; - } - - let last_ele_signum = last_ele.unwrap().y.signum(); - let ele_signum = ele.y.signum(); - - if last_ele_signum.is_nan() | ele_signum.is_nan() { - continue; - } - - if last_ele_signum != ele_signum { - // Do 50 iterations of newton's method, should be more than accurate - let x = { - let mut x1: f64 = last_ele.unwrap().x; - let mut x2: f64; - let mut fail: bool = false; - loop { - x2 = x1 - (self.function.get(x1) / self.function.get_derivative_1(x1)); - if !(self.min_x..self.max_x).contains(&x2) { - fail = true; - break; - } - - if (x2 - x1).abs() < resolution { - break; - } - - x1 = x2; - } - - match fail { - true => f64::NAN, - false => x1, - } - }; - - if !x.is_nan() { - root_list.push(Value::new(x, self.function.get(x))); - } - } - last_ele = Some(*ele); - } - self.output.roots = Some(root_list); + self.output.roots = Some( + newtons_method( + resolution, + self.min_x..self.max_x, + self.output.back.to_owned().unwrap(), + &|x: f64| self.function.get(x), + &|x: f64| self.function.get_derivative_1(x), + ) + .iter() + .map(|x| Value::new(*x, self.function.get(*x))) + .collect(), + ); } // Finds extrema fn extrema(&mut self) { let resolution: f64 = (self.pixel_width as f64 / (self.max_x - self.min_x).abs()) as f64; - let mut extrama_list: Vec = Vec::new(); - let mut last_ele: Option = None; - for ele in self.output.derivative.as_ref().unwrap().iter() { - if last_ele.is_none() { - last_ele = Some(*ele); - continue; - } - - let last_ele_signum = last_ele.unwrap().y.signum(); - let ele_signum = ele.y.signum(); - - if last_ele_signum.is_nan() | ele_signum.is_nan() { - continue; - } - - if last_ele_signum != ele_signum { - // Do 50 iterations of newton's method, should be more than accurate - let x = { - let mut x1: f64 = last_ele.unwrap().x; - let mut x2: f64; - let mut fail: bool = false; - loop { - x2 = x1 - - (self.function.get_derivative_1(x1) - / self.function.get_derivative_2(x1)); - if !(self.min_x..self.max_x).contains(&x2) { - fail = true; - break; - } - - if (x2 - x1).abs() < resolution { - break; - } - - x1 = x2; - } - - match fail { - true => f64::NAN, - false => x1, - } - }; - - if !x.is_nan() { - extrama_list.push(Value::new(x, self.function.get(x))); - } - } - last_ele = Some(*ele); - } - self.output.extrema = Some(extrama_list); + self.output.extrema = Some( + newtons_method( + resolution, + self.min_x..self.max_x, + self.output.derivative.to_owned().unwrap(), + &|x: f64| self.function.get_derivative_1(x), + &|x: f64| self.function.get_derivative_2(x), + ) + .iter() + .map(|x| Value::new(*x, self.function.get(*x))) + .collect(), + ); } pub fn display( diff --git a/src/misc.rs b/src/misc.rs index bb3cd1d..dc19f71 100644 --- a/src/misc.rs +++ b/src/misc.rs @@ -1,3 +1,7 @@ +use std::ops::Range; + +use eframe::egui::plot::Value; + cfg_if::cfg_if! { if #[cfg(target_arch = "wasm32")] { use wasm_bindgen::prelude::*; @@ -90,3 +94,66 @@ pub fn digits_precision(x: f64, digits: usize) -> f64 { let large_number: f64 = 10.0_f64.powf(digits as f64); (x * large_number).round() / large_number } + + +/// Implements newton's method of finding roots. +/// `threshold` is the target accuracy threshold +/// `range` is the range of valid x values (used to stop calculation when the point won't display anyways) +/// `data` is the data to iterate over (a Vector of egui's `Value` struct) +/// `f` is f(x) +/// `f_` is f'(x) +/// The function returns a list of `x` values where roots occur +pub fn newtons_method( + threshold: f64, range: Range, data: Vec, f: &dyn Fn(f64) -> f64, + f_1: &dyn Fn(f64) -> f64, +) -> Vec { + let mut output_list: Vec = Vec::new(); + let mut last_ele: Option = None; + for ele in data.iter() { + if last_ele.is_none() { + last_ele = Some(*ele); + continue; + } + + let last_ele_signum = last_ele.unwrap().y.signum(); + let ele_signum = ele.y.signum(); + + // If either are NaN, just continue iterating + if last_ele_signum.is_nan() | ele_signum.is_nan() { + continue; + } + + if last_ele_signum != ele_signum { + // Do 50 iterations of newton's method, should be more than accurate + let x = { + let mut x1: f64 = last_ele.unwrap().x; + let mut x2: f64; + let mut fail: bool = false; + loop { + x2 = x1 - (f(x1) / f_1(x1)); + if !range.contains(&x2) { + fail = true; + break; + } + + if (x2 - x1).abs() < threshold { + break; + } + + x1 = x2; + } + + match fail { + true => f64::NAN, + false => x1, + } + }; + + if !x.is_nan() { + output_list.push(x); + } + } + last_ele = Some(*ele); + } + output_list +}