linebender / kurbo

A Rust library for manipulating curves
Apache License 2.0
718 stars 69 forks source link

Possibility of using the `Float` trait instead of `f64`? #178

Open ctrlcctrlv opened 3 years ago

ctrlcctrlv commented 3 years ago

This relates to linebender/norad#108.

The num_traits crate declares a trait called Float which implements Neg, Add, Mul, Rem, Div, etc etc. It then implements this trait for f32 and f64.

pub trait Float: Num + Copy + NumCast + PartialOrd + Neg<Output = Self> {...}

Many libraries which work on numbers use this crate to be generic over floating point types.

I wonder if Kurbo could too?

Edited to add: I think this issue is fundamentally different from, but similar to, #159, because it basically lets the API consumer decide whether or not they are willing to take the risks of f32 lower precision.

raphlinus commented 3 years ago

My response is basically the same as #159. Generics make the types messier and potentially compile times longer and code more bloated (if there are two versions of the same algorithm). I think when interoperating with other code that prefers f32 (as is the case for Direct2D for example), putting in explicit precision-losing conversions is acceptable.

I do see that f64 is more precision than needed for font coordinates. There are other cases, though (scrolling is one) where the lower precision causes problems.

If I saw evidence of a real performance difference (as in #159), or evidence that ergonomics would be substantially improved, I would be more inclined to consider it. But in the absence of solid evidence, I am not really in favor.

cmyr commented 3 years ago

I remember we did some preliminary benchmarking a few years ago and found there was very little performance advantage to using f64, at least on x86_64. There would of course be some binary size advantages.

From where I stand right now, I have a good sense of the downsides of this change, but I don't really have a good sense of the upsides: I mean this in an honest way, I just haven't personally encountered any situations where I've really wished we were using f32, although I don't doubt that they exist.

my main concern with this is really the messier types.

There's a related argument that kurbo types should include an explicit generic parameter to represent an associated coordinate space; in this way you might have Point<PixelSpace> and Point<DesingSpace>, and you would have to explicitly move between them; this is a common thing to want in, say, a graphics editor. So far, though, this has also not felt worth the hassle.

Just to try and think this through, though, if we were going to do one or both of these, what I might imagine is having separate 'Raw' types that have the generic params, and then export simple types that have those params filled in, like:

pub struct RawPoint<T, R> {
    _space: PhantomData<R>,
    pub x: T,
    pub y: T,
}

pub trait CoordinateSpace: Copy + Eq {}

// everything is impl'd in terms of these params:
impl<T: CoordinateSpace, R: Float> std::ops::Sub<Point<T, R>> for Point<T, R> {
    type Output = Vec2<T, R>;
    #[inline]
    fn sub(self, other: Point<T, R>) -> Vec2<T, R> {
        Vec2::new(self.x - other.x, self.y - other.y)
    }
}

and then in lib.rs, do something like:

pub type Point = RawPoint<GenericCoordinateSpace, f64>;

And then hopefully the user-facing code looks basically the same as it does now, unless you happen to need to customize.

My big problem with this is that I think the docs would end up really gross, and I don't have any idea for how we might address that...

JAicewizard commented 3 years ago

I think that maybe when posits start to be implemented in consumer CPUs there would be a bigger argument for making this generic. However, we are looking 3/4+ years into the future. I myself dont see the point. I just had to debug performance issues related to num_trait, so when not being careful there certainly are also downsides to making this generic.

richard-uk1 commented 1 year ago

I'm very much against making f64 generic, because I'm currently involved in a project that uses the euclid crate and the compile times are painful.

ctrlcctrlv commented 1 year ago

Couldn't it be optional? Please see integer_or_float for how I handle that.

ctrlcctrlv commented 1 year ago

https://docs.rs/integer_or_float/latest/integer_or_float/backing_types/constant.IOF_BACKING_STORE_BITLEN.html https://docs.rs/integer_or_float/latest/integer_or_float/backing_types/type.f_iof.html

It'd be no problem in my understanding to do something like this in kurbo?

ctrlcctrlv commented 1 year ago

Sorry to triple post. Kept thinking of more parallels between the problems.

IOF has a lot of feature flags. One is num-traits.

That enables the entire tree under mod num_traits_impl:

https://docs.rs/crate/integer_or_float/latest/source/src/num_traits_impl/float.rs

I use it in MFEK but it can be turned off.

ctrlcctrlv commented 1 year ago

@raphlinus I understand IOF has been rejected for kurbo and norad, no problem. But if it'd help us come to consensus, I can split off the x64_backing_store feature to a new crate to be called MFEK/generic-primitive.rlib. That crate will then be the one to have num_traits_impl feature. Does this work for @/everyone?

cmyr commented 1 year ago

So I think using a feature to alternate between 32-and-64 bit numbers isn't really an option in a library crate; features are supposed to be additive, not to toggle between alternate incompatible behaviours. This would be a concrete problem: if we imagine a float64 feature that swaps in f64 for f32, if I had a crate without this feature that was assuming kurbo always used f32, and then it had a dependency that used kurbo with this feature, our crate would no longer compile.

For the latter question, It isn't clear to me what the proposal is? What are we trying to reach consensus on?

ctrlcctrlv commented 1 year ago

if I had a crate without this feature that was assuming kurbo always used f32

That'd be a crate that is still working only because it has not yet updated to use the kurbo::float type which would come from generic_primitives::float compiled with x64_backing_store.

richard-uk1 commented 1 year ago

I think the neatest solution is to have a kurbo_f32 crate for the 32-bit case. In the future it could potentially be added into kurbo proper as the module kurbo::f32 and behind a feature flag. Would this solve your issues @ctrlcctrlv?

You could create it by cloning the kurbo crate and running a replace f64 -> f32, and then fix up anything manually that you need to.

ctrlcctrlv commented 1 year ago

@derekdreery It would.

richard-uk1 commented 1 year ago

We talked about this at the office hours. The consensus was that on CPU, there really isn't much perf benefit to using f32 over f64, so the suggested workflow for people wanting f32 is to convert up when calling kurbo fns, and then cast down (as f32) the output.

This doesn't stop anyone from making an f32 version of kurbo out-of-tree.

ctrlcctrlv commented 1 year ago

Please reopen this because I'm trying and you're not correct. Tests fail, one loops forever.

ctrlcctrlv commented 1 year ago

Possibly related to #50.

richard-uk1 commented 1 year ago

Please reopen this because I'm trying and you're not correct. Tests fail, one loops forever.

Can you say more? What are you trying?

ctrlcctrlv commented 1 year ago

The problem is certainly related to choice of epsila.

diff --git a/src/common.rs b/src/common.rs
index f4e9107..07944b4 100644
--- a/src/common.rs
+++ b/src/common.rs
@@ -286,7 +286,7 @@ pub fn solve_itp(
 ) -> f32 {
     let n1_2 = (((b - a) / epsilon).log2().ceil() - 1.0).max(0.0) as usize;
     let nmax = n0 + n1_2;
-    let mut scaled_epsilon = epsilon * (1u64 << nmax) as f32;
+    let mut scaled_epsilon = epsilon * (1u32 << nmax) as f32;
     while b - a > 2.0 * epsilon {
         let x1_2 = 0.5 * (a + b);
         let r = scaled_epsilon - 0.5 * (b - a);
@@ -540,7 +540,7 @@ mod tests {

     fn verify<const N: usize>(mut roots: ArrayVec<f32, N>, expected: &[f32]) {
         assert_eq!(expected.len(), roots.len());
-        let epsilon = 1e-12;
+        let epsilon = f32::EPSILON;
         roots.sort_by(|a, b| a.partial_cmp(b).unwrap());
         for i in 0..expected.len() {
             assert!((roots[i] - expected[i]).abs() < epsilon);
ctrlcctrlv commented 1 year ago

with that patch my CPU thrashing stopped and now i just get failures

fred@デブ狸~/Workspace/kurbo% cargo test
    Finished test [unoptimized + debuginfo] target(s) in 0.01s
     Running unittests src/lib.rs (target/debug/deps/kurbo-38cdea1756b97c3e)

running 75 tests
test affine::tests::affine_inv ... FAILED
test affine::tests::affine_mul ... FAILED
test bezpath::tests::test_contains ... ok
test affine::tests::affine_basic ... FAILED
test bezpath::tests::test_elements_to_segments_closepath_refers_to_last_moveto ... ok
test bezpath::tests::test_elements_to_segments_starts_on_closepath - should panic ... ok
test bezpath::tests::test_get_seg ... ok
test bezpath::tests::test_intersect_cubic ... FAILED
test bezpath::tests::test_intersect_line ... ok
test bezpath::tests::test_intersect_qad ... ok
test bezpath::tests::test_must_not_start_on_quad - should panic ... ok
test common::tests::test_inv_arclen ... ok
test circle::tests::area_sign ... FAILED
test common::tests::test_solve_cubic ... FAILED
test common::tests::test_solve_itp ... FAILED
test common::tests::test_solve_quadratic ... ok
test cubicbez::tests::cubicbez_arclen ... ok
test cubicbez::tests::cubicbez_approx_spline ... FAILED
test cubicbez::tests::cubicbez_cubics_to_quadratic_splines ... ok
test cubicbez::tests::cubicbez_deriv ... FAILED
test cubicbez::tests::cubicbez_extrema ... ok
test cubicbez::tests::cubicbez_inflections ... ok
test cubicbez::tests::cubicbez_inv_arclen_accuracy ... FAILED
test cubicbez::tests::cubicbez_signed_area ... FAILED
test cubicbez::tests::cubicbez_signed_area_linear ... FAILED
test cubicbez::tests::degenerate_to_quads ... ok
test line::tests::line_arclen ... ok
test ellipse::tests::area_sign ... FAILED
test line::tests::line_is_finite ... ok
test cubicbez::tests::cubicbez_toquads ... FAILED
test mindist::tests::test_choose ... ok
test mindist::tests::test_d_rk ... ok
test cubicbez::tests::cubicbez_nearest ... FAILED
test point::tests::display ... ok
test point::tests::distance ... ok
test quadbez::tests::quadbez_arclen ... ok
test point::tests::point_arithmetic ... ok
test quadbez::tests::quadbez_deriv ... FAILED
test quadbez::tests::quadbez_arclen_pathological ... FAILED
test quadbez::tests::quadbez_extrema ... ok
test mindist::tests::test_overflow ... ok
test quadbez::tests::quadbez_nearest ... ok
test quadbez::tests::quadbez_raise ... FAILED
test quadbez::tests::quadbez_nearest_low_order ... ok
test quadbez::tests::quadbez_subsegment ... FAILED
test quadbez::tests::quadbez_signed_area ... FAILED
test quadspline::tests::four_points_implicit_on_curve ... ok
test quadspline::tests::no_points_no_quads ... ok
test quadspline::tests::one_point_no_quads ... ok
test cubicbez::tests::cubicbez_inv_arclen ... FAILED
test quadspline::tests::three_points_same_quad ... ok
test rect::tests::aspect_ratio ... ok
test quadspline::tests::two_points_no_quads ... ok
test rect::tests::area_sign ... ok
test rect::tests::contained_rect_with_aspect_ratio ... ok
test rect::tests::display ... FAILED
test rounded_rect::tests::area ... FAILED
test rounded_rect::tests::winding ... ok
test rounded_rect::tests::bez_conversion ... FAILED
test size::tests::aspect_ratio ... ok
test size::tests::display ... ok
test svg::tests::test_parse_svg ... ok
test svg::tests::test_parse_svg2 ... ok
test svg::tests::test_parse_svg_arc ... FAILED
test svg::tests::test_parse_svg_arc_pie ... ok
test svg::tests::test_write_svg_single ... ok
test svg::tests::test_write_svg_two_move ... ok
test translate_scale::tests::conversions ... ok
test svg::tests::test_write_svg_two_nomove ... ok
test translate_scale::tests::inverse ... ok
test translate_scale::tests::translate_scale ... ok
test vec2::tests::display ... ok
test mindist::tests::test_out_of_order ... ok
test svg::tests::test_serialize_deserialize ... ok
test mindist::tests::test_mindist ... ok

failures:

---- affine::tests::affine_inv stdout ----
thread 'affine::tests::affine_inv' panicked at '(0.0, 1.0000005) != (0.0, 1.0)', src/affine.rs:333:9

---- affine::tests::affine_mul stdout ----
thread 'affine::tests::affine_mul' panicked at '(30.000002, 42.4) != (30.0, 42.4)', src/affine.rs:333:9
note: run with `RUST_BACKTRACE=1` environment variable to display a backtrace

---- affine::tests::affine_basic stdout ----
thread 'affine::tests::affine_basic' panicked at '(-4.0, 2.9999998) != (-4.0, 3.0)', src/affine.rs:333:9

---- bezpath::tests::test_intersect_cubic stdout ----
thread 'bezpath::tests::test_intersect_cubic' panicked at '0.59259254 != 0.5925926', src/bezpath.rs:1341:9

---- circle::tests::area_sign stdout ----
thread 'circle::tests::area_sign' panicked at '78.53982 != 78.53984', src/circle.rs:366:9

---- common::tests::test_solve_cubic stdout ----
thread 'common::tests::test_solve_cubic' panicked at 'assertion failed: (roots[i] - expected[i]).abs() < epsilon', src/common.rs:546:13

---- common::tests::test_solve_itp stdout ----
thread 'common::tests::test_solve_itp' panicked at 'attempt to shift left with overflow', src/common.rs:289:40

---- cubicbez::tests::cubicbez_approx_spline stdout ----
thread 'cubicbez::tests::cubicbez_approx_spline' panicked at 'assertion failed: `(left == right)`
  left: `QuadBez { p0: (550.0, 258.0), p1: (1673.6658, 767.5165), p2: (1934.0, 1554.0) }`,
 right: `QuadBez { p0: (550.0, 258.0), p1: (1673.6658, 767.5164), p2: (1934.0, 1554.0) }`', src/cubicbez.rs:848:9

---- cubicbez::tests::cubicbez_deriv stdout ----
thread 'cubicbez::tests::cubicbez_deriv' panicked at 'assertion failed: (d - d_approx).hypot() < delta * 2.0', src/cubicbez.rs:668:13

---- cubicbez::tests::cubicbez_inv_arclen_accuracy stdout ----
thread 'cubicbez::tests::cubicbez_inv_arclen_accuracy' panicked at 'attempt to shift left with overflow', src/common.rs:289:40

---- cubicbez::tests::cubicbez_signed_area stdout ----
thread 'cubicbez::tests::cubicbez_signed_area' panicked at 'assertion failed: ((Affine::rotate(0.5) * c).signed_area() - 0.75).abs() < epsilon', src/cubicbez.rs:757:9

---- cubicbez::tests::cubicbez_signed_area_linear stdout ----
thread 'cubicbez::tests::cubicbez_signed_area_linear' panicked at 'assertion failed: ((Affine::translate((1.0, 0.0)) * c).signed_area() - 1.0).abs() < epsilon', src/cubicbez.rs:748:9

---- ellipse::tests::area_sign stdout ----
thread 'ellipse::tests::area_sign' panicked at '78.5398 != 78.53982', src/ellipse.rs:282:9

---- cubicbez::tests::cubicbez_toquads stdout ----
thread 'cubicbez::tests::cubicbez_toquads' panicked at 'got 0.000001013279 wanted 0.0000010000002', src/cubicbez.rs:826:21

---- cubicbez::tests::cubicbez_nearest stdout ----
thread 'cubicbez::tests::cubicbez_nearest' panicked at 'got Nearest { distance_sq: 5.47562e-12, t: 0.6999984 } expected 0.7', src/cubicbez.rs:765:13

---- quadbez::tests::quadbez_deriv stdout ----
thread 'quadbez::tests::quadbez_deriv' panicked at 'assertion failed: (d - d_approx).hypot() < delta * 2.0', src/quadbez.rs:410:13

---- quadbez::tests::quadbez_arclen_pathological stdout ----
thread 'quadbez::tests::quadbez_arclen_pathological' panicked at '2.0008736 != 2.0008738', src/quadbez.rs:432:9

---- quadbez::tests::quadbez_raise stdout ----
thread 'quadbez::tests::quadbez_raise' panicked at '(3.6259997, 3.8469996) != (3.6259995, 3.8469994)', src/quadbez.rs:394:9

---- quadbez::tests::quadbez_subsegment stdout ----
thread 'quadbez::tests::quadbez_subsegment' panicked at '(3.9537396, 3.7258298) != (3.9537396, 3.7258296)', src/quadbez.rs:394:9

---- quadbez::tests::quadbez_signed_area stdout ----
thread 'quadbez::tests::quadbez_signed_area' panicked at 'assertion failed: ((Affine::rotate(0.5) * q).signed_area() - 2.0 / 3.0).abs() < epsilon', src/quadbez.rs:474:9

---- cubicbez::tests::cubicbez_inv_arclen stdout ----
thread 'cubicbez::tests::cubicbez_inv_arclen' panicked at 'at accuracy 1.0000002e-5, wanted 103.52598 got 103.526', src/cubicbez.rs:706:17

---- rect::tests::display stdout ----
thread 'rect::tests::display' panicked at 'assertion failed: `(left == right)`
  left: `"Rect { (10, 12.23214) (22.222221×23.099998) }"`,
 right: `"Rect { (10, 12.23214) (22.222222222×23.1) }"`', src/rect.rs:716:9

---- rounded_rect::tests::area stdout ----
thread 'rounded_rect::tests::area' panicked at 'assertion failed: (circle.area() - rounded_rect.area()).abs() < epsilon', src/rounded_rect.rs:451:9

---- rounded_rect::tests::bez_conversion stdout ----
thread 'rounded_rect::tests::bez_conversion' panicked at 'assertion failed: (rect.area() - p.area()).abs() < epsilon', src/rounded_rect.rs:474:9

---- svg::tests::test_parse_svg_arc stdout ----
thread 'svg::tests::test_parse_svg_arc' panicked at 'assertion failed: `(left == right)`
  left: `4`,
 right: `3`', src/svg.rs:507:9

failures:
    affine::tests::affine_basic
    affine::tests::affine_inv
    affine::tests::affine_mul
    bezpath::tests::test_intersect_cubic
    circle::tests::area_sign
    common::tests::test_solve_cubic
    common::tests::test_solve_itp
    cubicbez::tests::cubicbez_approx_spline
    cubicbez::tests::cubicbez_deriv
    cubicbez::tests::cubicbez_inv_arclen
    cubicbez::tests::cubicbez_inv_arclen_accuracy
    cubicbez::tests::cubicbez_nearest
    cubicbez::tests::cubicbez_signed_area
    cubicbez::tests::cubicbez_signed_area_linear
    cubicbez::tests::cubicbez_toquads
    ellipse::tests::area_sign
    quadbez::tests::quadbez_arclen_pathological
    quadbez::tests::quadbez_deriv
    quadbez::tests::quadbez_raise
    quadbez::tests::quadbez_signed_area
    quadbez::tests::quadbez_subsegment
    rect::tests::display
    rounded_rect::tests::area
    rounded_rect::tests::bez_conversion
    svg::tests::test_parse_svg_arc

test result: FAILED. 50 passed; 25 failed; 0 ignored; 0 measured; 0 filtered out; finished in 0.22s
ctrlcctrlv commented 1 year ago

All epsila need to be calculated using f32::EPSILON or f64::EPSILON as their basis, all these lines need changed:

fred@デブ狸~/Workspace/kurbo% rg 'epsilon[ \t]*='
src/line.rs
301:        let epsilon = 1e-9;

// Note: I already changed these in my version but they're wrong in master.
src/common.rs
289:    let mut scaled_epsilon = epsilon * (1u32 << nmax) as f32;
543:        let epsilon = f32::EPSILON;

src/rounded_rect_radii.rs
96:        let epsilon = 1e-9;

src/rounded_rect.rs
441:        let epsilon = 1e-9;
473:        let epsilon = 1e-7;

src/quadbez.rs
444:        let epsilon = 1e-12;
459:        let epsilon = 1e-12;
472:        let epsilon = 1e-12;

src/cubicbez.rs
744:        let epsilon = 1e-12;
755:        let epsilon = 1e-12;
817:                let epsilon = 1e-12;
ctrlcctrlv commented 1 year ago

If we need to raise them by powers we should multiply fxx::EPSILON by (10.0).powf(n) I think.

ctrlcctrlv commented 1 year ago

Here's a diff that fixes some affine tests

diff --git a/src/affine.rs b/src/affine.rs
index 00d770e..d403e76 100644
--- a/src/affine.rs
+++ b/src/affine.rs
@@ -330,7 +330,7 @@ mod tests {
     use std::f32::consts::PI;

     fn assert_near(p0: Point, p1: Point) {
-        assert!((p1 - p0).hypot() < 1e-9, "{p0:?} != {p1:?}");
+        assert!((p1 - p0).hypot() < f32::EPSILON * (10f32.powi(3)), "{p0:?} != {p1:?}");
     }

     #[test]
ctrlcctrlv commented 1 year ago

I added this as src/epsilon.rs

use core::ops::Mul;

const EPSILON: f64 = f64::EPSILON;
const MAX_10_EXP: i32 = f64::MAX_10_EXP;

/// `kurbo`'s epsilon (ε).
#[repr(packed)]
#[derive(Copy, Clone, Debug)]
pub struct Epsilon {
    /// Value of this epsilon. By default, [`f64::EPSILON`].
    pub value: f64,
    /// Magnitude of this epsilon. By default, [`f64::MAX_10_EXP`].
    pub magnitude: i32,
}

impl Default for Epsilon {
    #[inline]
    fn default() -> Epsilon {
        Self::new()
    }
}

impl Epsilon {
    /// Raises or lowers the magnitude of `kurbo`'s epsilon by `magnitude`, returning a new
    /// [`Epsilon`].
    ///
    /// Returns a new `Epsilon` representing **ε × 10_ᵐ_** where _m_ = magnitude.
    #[inline]
    pub const fn ten_pow(magnitude: i32) -> Self {
        debug_assert!(MAX_10_EXP + magnitude <= MAX_10_EXP);

        Epsilon {
            magnitude: MAX_10_EXP + magnitude,
            value: f64::mul(EPSILON, 10.0 * magnitude as f64)
        }
    }
}

impl Epsilon {
    /// Create `kurbo`'s default epsilon.
    #[inline]
    pub const fn new() -> Self {
        Epsilon {
            value: EPSILON,
            magnitude: MAX_10_EXP
        }
    }
}
ctrlcctrlv commented 1 year ago

I'd also like kurbo to start declaring a pub const, pub const FLOAT_WIDTH: u16 = 64.

In my version we'll declare that as 32.

That's because I will tell people in docs to use kurbo_f32 as kurbo for ease of use...

richard-uk1 commented 1 year ago

@ctrlcctrlv what would the purpose of the pub const be?

ctrlcctrlv commented 1 year ago

Knowing which kurbo I'm using without going through tangled mess of cargo tree.

richard-uk1 commented 1 year ago

We talked about this in the office hours again. @raphlinus was concerned that some of the stability analysis of the algorithms would become more important for f32, and that this would require a lot of work. This is in addition to the hit we'd take on compile times. We need to understand what the motivation is for having 32bit versions of kurbo objects. Is it

For example, if there were use cases like wanting cubic Bézier curves in memory for the GPU, we could add a CubicBez32 struct, where all the algorithms still run in f64 and the result is cast at the end.

edit tagged Raph because I'm speaking on his behalf and want to check that I'm faithfully representing him.

richard-uk1 commented 1 year ago

One possible use of f32 on the CPU side that came up at my work is in R-trees, for keeping the size of the index small.

ctrlcctrlv commented 1 year ago

For example, if there were use cases like wanting cubic Bézier curves in memory for the GPU, we could add a CubicBez32 struct, where all the algorithms still run in f64 and the result is cast at the end.

Actually, I find this solution perfectly reasonable. Do that.