From 8b62a2bd8c1e44db401ed4ec04d4f9fbf7674e63 Mon Sep 17 00:00:00 2001 From: Michael Kirk Date: Fri, 26 Jan 2024 16:50:12 -0800 Subject: [PATCH] use libm remquo and copysign for a little speedup MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit This makes our implement of sincosd look more like the cpp version. $ cargo bench --bench="*" -- --baseline=4891fc4 Finished bench [optimized] target(s) in 0.02s Running benches/geodesic_benchmark.rs (target/release/deps/geodesic_benchmark-cda172d0984b8505) direct (c wrapper)/default time: [24.046 µs 24.071 µs 24.099 µs] change: [-0.2131% +0.0210% +0.2593%] (p = 0.86 > 0.05) No change in performance detected. Found 9 outliers among 100 measurements (9.00%) 3 (3.00%) low severe 2 (2.00%) low mild 1 (1.00%) high mild 3 (3.00%) high severe direct (rust impl)/default time: [26.129 µs 26.168 µs 26.211 µs] change: [-5.5845% -5.3729% -5.1792%] (p = 0.00 < 0.05) Performance has improved. Found 9 outliers among 100 measurements (9.00%) 4 (4.00%) high mild 5 (5.00%) high severe inverse (c wrapper)/default time: [45.061 µs 45.141 µs 45.227 µs] change: [+0.3738% +0.6326% +0.8919%] (p = 0.00 < 0.05) Change within noise threshold. Found 7 outliers among 100 measurements (7.00%) 3 (3.00%) high mild 4 (4.00%) high severe inverse (rust impl)/default time: [67.739 µs 67.796 µs 67.865 µs] change: [-3.2442% -3.0580% -2.8479%] (p = 0.00 < 0.05) Performance has improved. Found 11 outliers among 100 measurements (11.00%) 1 (1.00%) low severe 5 (5.00%) low mild 2 (2.00%) high mild 3 (3.00%) high severe --- Cargo.toml | 1 + README.md | 10 ++++---- src/geomath.rs | 64 ++++++++++---------------------------------------- 3 files changed, 18 insertions(+), 57 deletions(-) diff --git a/Cargo.toml b/Cargo.toml index e765674..375eb17 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -25,6 +25,7 @@ default = ["accurate"] [dependencies] accurate = { version = "0.3", optional = true, default-features = false } +libm = { version = "0.2.8", default-features = false } [dev-dependencies] approx = "0.5.1" diff --git a/README.md b/README.md index bdc9019..68f636e 100644 --- a/README.md +++ b/README.md @@ -87,16 +87,16 @@ Which produces output like: ```text direct (c wrapper)/default - time: [24.055 µs 24.085 µs 24.117 µs] + time: [24.046 µs 24.071 µs 24.099 µs] direct (rust impl)/default - time: [27.760 µs 27.810 µs 27.867 µs] + time: [26.129 µs 26.168 µs 26.211 µs] inverse (c wrapper)/default - time: [46.461 µs 47.435 µs 48.557 µs] + time: [45.061 µs 45.141 µs 45.227 µs] inverse (rust impl)/default - time: [70.488 µs 70.841 µs 71.356 µs] + time: [67.739 µs 67.796 µs 67.865 µs] ``` -Showing that, at least in this benchmark, the Rust implementation is 16-52% slower than the c bindings. +Showing that, at least in this benchmark, the Rust implementation is 10-50% slower than the c bindings. diff --git a/src/geomath.rs b/src/geomath.rs index 12bf174..4482f06 100644 --- a/src/geomath.rs +++ b/src/geomath.rs @@ -136,66 +136,26 @@ pub fn ang_diff(x: f64, y: f64) -> (f64, f64) { } } -pub fn fmod(x: f64, y: f64) -> f64 { - x % y -} - /// Compute sine and cosine of x in degrees pub fn sincosd(x: f64) -> (f64, f64) { - // r = math.fmod(x, 360) if Math.isfinite(x) else Math.nan - let mut r = if x.is_finite() { - fmod(x, 360.0) - } else { - std::f64::NAN - }; - - // q = 0 if Math.isnan(r) else int(round(r / 90)) - let mut q = if r.is_nan() { - 0 - } else { - (r / 90.0).round() as i32 - }; + let (mut r, q) = libm::remquo(x, 90.0); - // r -= 90 * q; r = math.radians(r) - r -= 90.0 * q as f64; r = r.to_radians(); - // s = math.sin(r); c = math.cos(r) - let s = r.sin(); - let c = r.cos(); + let (mut sinx, mut cosx) = r.sin_cos(); - // q = q % 4 - q %= 4; - - // if q == 1: - // s, c = c, -s - // elif q == 2: - // s, c = -s, -c - // elif q == 3: - // s, c = -c, s - - let q = if q < 0 { q + 4 } else { q }; - - let (s, c) = if q == 1 { - (c, -s) - } else if q == 2 { - (-s, -c) - } else if q == 3 { - (-c, s) - } else { - debug_assert_eq!(q, 0); - (s, c) + (sinx, cosx) = match q as u32 & 3 { + 0 => (sinx, cosx), + 1 => (cosx, -sinx), + 2 => (-sinx, -cosx), + 3 => (-cosx, sinx), + _ => unreachable!(), }; - // # Remove the minus sign on -0.0 except for sin(-0.0). - // # On Windows 32-bit with python 2.7, math.fmod(-0.0, 360) = +0.0 - // # (x, c) here fixes this bug. See also Math::sincosd in the C++ library. - // # AngNormalize has a similar fix. - // s, c = (x, c) if x == 0 else (0.0+s, 0.0+c) - // return s, c - let (s, c) = if x == 0.0 { (x, c) } else { (0.0 + s, 0.0 + c) }; - - (s, c) + if sinx == 0.0 { + sinx = libm::copysign(sinx, x); // special values from F.10.1.13 + } + (sinx, cosx) } // Compute atan2(y, x) with result in degrees