From e1f3d3cd80fdd6e214fbdcc14e9afeb6633abab5 Mon Sep 17 00:00:00 2001 From: Tinko Bartels Date: Tue, 23 May 2023 20:29:23 +0200 Subject: [PATCH] Use error bound from Ozaki et al in orient2d, remove (now) unused constant CCWERRBOUND_A. --- src/lib.rs | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/src/lib.rs b/src/lib.rs index 85be680..28137e0 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -51,7 +51,6 @@ pub struct Coord3D> { const SPLITTER: f64 = 134_217_729f64; const EPSILON: f64 = 0.000_000_000_000_000_111_022_302_462_515_65; const RESULTERRBOUND: f64 = (3.0 + 8.0 * EPSILON) * EPSILON; -const CCWERRBOUND_A: f64 = (3.0 + 16.0 * EPSILON) * EPSILON; const CCWERRBOUND_B: f64 = (2.0 + 12.0 * EPSILON) * EPSILON; const CCWERRBOUND_C: f64 = (9.0 + 64.0 * EPSILON) * EPSILON * EPSILON; const O3DERRBOUND_A: f64 = (7.0 + 56.0 * EPSILON) * EPSILON; @@ -86,11 +85,14 @@ pub fn orient2d>(pa: Coord, pb: Coord, pc: Coord) -> f64 { let detright = (pa.y - pc.y) * (pb.x - pc.x); let det = detleft - detright; - // The detsum calculation was changed to require only one branch on the likely execution path. + // The errbound calculation was changed to require only one branch on the likely execution path. // This improves performance on modern processors as discussed by Ozaki et al in // https://doi.org/10.1007/s10543-015-0574-9 - let detsum = abs(detleft) + abs(detright); - let errbound = CCWERRBOUND_A * detsum; + // The underflow guard "+ u_N" was omitted because orient2dadapt(...) would not guarantee + // correct results in cases of underflow, the derivation of THETA is given in the reference. + let detsum = abs(detleft + detright); + const THETA: f64 = 3.3306690621773722e-16; + let errbound = THETA * detsum; if det >= errbound || -det >= errbound { det } else {