/rust/registry/src/index.crates.io-1949cf8c6b5b557f/polycool-0.4.0/src/yuksel.rs
Line | Count | Source |
1 | | // Copyright 2025 the Kurbo Authors |
2 | | // SPDX-License-Identifier: Apache-2.0 OR MIT |
3 | | |
4 | | //! An implementation of Yuksel's robust version of Newton's algorithm. |
5 | | |
6 | 0 | fn different_signs(x: f64, y: f64) -> bool { |
7 | 0 | (x < 0.0) != (y < 0.0) |
8 | 0 | } |
9 | | |
10 | 0 | pub(crate) fn find_root<F: Fn(f64) -> f64, DF: Fn(f64) -> f64>( |
11 | 0 | f: F, |
12 | 0 | deriv: DF, |
13 | 0 | mut lower: f64, |
14 | 0 | mut upper: f64, |
15 | 0 | val_lower: f64, |
16 | 0 | val_upper: f64, |
17 | 0 | x_error: f64, |
18 | 0 | ) -> f64 { |
19 | 0 | if !val_lower.is_finite() || !val_upper.is_finite() { |
20 | 0 | return f64::NAN; |
21 | 0 | } |
22 | 0 | debug_assert!(different_signs(val_lower, val_upper)); |
23 | | |
24 | 0 | let mut x = lower + (upper - lower) / 2.0; |
25 | 0 | let mut step = (upper - lower) / 2.0; |
26 | | |
27 | 0 | if step.abs() <= x_error { |
28 | 0 | return x; |
29 | 0 | } |
30 | | |
31 | 0 | while step.abs() > x_error && x.is_finite() { |
32 | | // There's potentially a performance gain to be had by evaluating the |
33 | | // derivative and the polynomial "jointly", so that subexpressions can |
34 | | // be shared. At least for the cubic case, the compiler is smart enough |
35 | | // to do that for us. |
36 | 0 | let deriv_x = deriv(x); |
37 | 0 | let val_x = f(x); |
38 | | |
39 | | // By returning early on zero, we make it impossible for `step` to |
40 | | // be NaN. This seems to benchmark a little better than checking for |
41 | | // NaN after the fact. |
42 | 0 | if val_x == 0.0 { |
43 | 0 | return x; |
44 | 0 | } |
45 | 0 | let root_in_first_half = different_signs(val_lower, val_x); |
46 | 0 | if root_in_first_half { |
47 | 0 | upper = x; |
48 | 0 | } else { |
49 | 0 | lower = x; |
50 | 0 | } |
51 | | |
52 | 0 | debug_assert!(deriv_x.is_finite()); |
53 | 0 | debug_assert!(val_x.is_finite()); |
54 | | |
55 | 0 | step = -val_x / deriv_x; |
56 | 0 | let mut new_x = x + step; |
57 | | |
58 | 0 | if new_x <= lower || new_x >= upper { |
59 | 0 | new_x = lower + (upper - lower) / 2.0; |
60 | | |
61 | 0 | if new_x == upper || new_x == lower { |
62 | | // This should be rare, but it happens if they ask for more |
63 | | // accuracy than is reasonable. For example, suppose (because |
64 | | // of large coefficients) the output value jumps from -1.0 |
65 | | // to 1.0 between adjacent floats and they ask for an output |
66 | | // error of smaller than 0.5. Then we'll eventually shrink |
67 | | // the search interval to a pair of adjacent floats and hit |
68 | | // this case. |
69 | 0 | return new_x; |
70 | 0 | } |
71 | 0 | } |
72 | 0 | step = new_x - x; |
73 | 0 | x = new_x; |
74 | | } |
75 | 0 | x |
76 | 0 | } Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<5>>::roots_between_with_buffer<5>::{closure#0}, <polycool::poly::Poly<5>>::roots_between_with_buffer<5>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<5>>::roots_between_with_buffer<4>::{closure#0}, <polycool::poly::Poly<5>>::roots_between_with_buffer<4>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<5>>::roots_between_with_buffer<6>::{closure#0}, <polycool::poly::Poly<5>>::roots_between_with_buffer<6>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<5>>::roots_between_with_buffer<7>::{closure#0}, <polycool::poly::Poly<5>>::roots_between_with_buffer<7>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<5>>::roots_between_with_buffer<8>::{closure#0}, <polycool::poly::Poly<5>>::roots_between_with_buffer<8>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<5>>::roots_between_with_buffer<9>::{closure#0}, <polycool::poly::Poly<5>>::roots_between_with_buffer<9>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<6>>::roots_between_with_buffer<6>::{closure#0}, <polycool::poly::Poly<6>>::roots_between_with_buffer<6>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<6>>::roots_between_with_buffer<5>::{closure#0}, <polycool::poly::Poly<6>>::roots_between_with_buffer<5>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<6>>::roots_between_with_buffer<7>::{closure#0}, <polycool::poly::Poly<6>>::roots_between_with_buffer<7>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<6>>::roots_between_with_buffer<8>::{closure#0}, <polycool::poly::Poly<6>>::roots_between_with_buffer<8>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<6>>::roots_between_with_buffer<9>::{closure#0}, <polycool::poly::Poly<6>>::roots_between_with_buffer<9>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<7>>::roots_between_with_buffer<7>::{closure#0}, <polycool::poly::Poly<7>>::roots_between_with_buffer<7>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<7>>::roots_between_with_buffer<6>::{closure#0}, <polycool::poly::Poly<7>>::roots_between_with_buffer<6>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<7>>::roots_between_with_buffer<8>::{closure#0}, <polycool::poly::Poly<7>>::roots_between_with_buffer<8>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<7>>::roots_between_with_buffer<9>::{closure#0}, <polycool::poly::Poly<7>>::roots_between_with_buffer<9>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<8>>::roots_between_with_buffer<8>::{closure#0}, <polycool::poly::Poly<8>>::roots_between_with_buffer<8>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<8>>::roots_between_with_buffer<7>::{closure#0}, <polycool::poly::Poly<8>>::roots_between_with_buffer<7>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<8>>::roots_between_with_buffer<9>::{closure#0}, <polycool::poly::Poly<8>>::roots_between_with_buffer<9>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<9>>::roots_between_with_buffer<9>::{closure#0}, <polycool::poly::Poly<9>>::roots_between_with_buffer<9>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<9>>::roots_between_with_buffer<8>::{closure#0}, <polycool::poly::Poly<9>>::roots_between_with_buffer<8>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<10>>::roots_between_with_buffer<9>::{closure#0}, <polycool::poly::Poly<10>>::roots_between_with_buffer<9>::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly::Poly<4>>::one_root::{closure#0}, <polycool::poly::Poly<4>>::one_root::{closure#1}>Unexecuted instantiation: polycool::yuksel::find_root::<<polycool::poly_dyn::PolyDyn>::roots_between_with_buffer::{closure#1}, <polycool::poly_dyn::PolyDyn>::roots_between_with_buffer::{closure#2}> |