zoe/math/
float.rs

1#![allow(dead_code)]
2
3use std::{
4    fmt::Debug,
5    ops::{Add, AddAssign, Div, DivAssign, Mul, Neg, Sub},
6    simd::SimdElement,
7    str::FromStr,
8};
9
10use crate::private::Sealed;
11
12/// Takes square root using the Babylonian Method using 40 iterations.
13// Relevant discussion: <https://www.codeproject.com/Articles/69941/Best-Square-Root-Method-Algorithm-Function-Precisi>
14pub(crate) const fn sqrt_baby(n: f64) -> f64 {
15    let mut i: f64 = 0.0;
16    while i * i <= n {
17        i += 0.1_f64;
18    }
19
20    let mut x1 = i;
21    let mut x2 = 0.0;
22
23    let mut j: usize = 0;
24    while j < 40 {
25        x2 = (n / x1).midpoint(x1);
26        x1 = x2;
27        j += 1;
28    }
29    x2
30}
31
32/// Trait for providing generic functionality over floating point numbers.
33pub trait Float:
34    Sub<Output = Self>
35    + Add<Output = Self>
36    + Mul<Output = Self>
37    + Div<Output = Self>
38    + Neg<Output = Self>
39    + AddAssign
40    + DivAssign
41    + Sized
42    + Default
43    + PartialEq
44    + PartialOrd
45    + Copy
46    + Debug
47    + SimdElement
48    + Sized
49    + FromStr
50    + std::iter::Sum<Self>
51    + for<'a> std::iter::Sum<&'a Self>
52    + Sealed {
53    const MIN: Self;
54    const MIN_POSITIVE: Self;
55    const MAX: Self;
56    const ZERO: Self;
57    const ONE: Self;
58    const INFINITY: Self;
59
60    /// Generic absolute value for [`Float`]
61    fn abs(self) -> Self;
62    /// Generic minimum of 2 values for [`Float`]
63    fn min(self, other: Self) -> Self;
64    /// Generic log2 for [`Float`]
65    fn log2(self) -> Self;
66    /// Generic `is_nan` calculation for [`Float`]
67    fn is_nan(self) -> bool;
68    /// Generic `is_infinite` for [`Float`]
69    fn is_infinite(self) -> bool;
70    /// Generic natural logarithm for [`Float`]
71    fn ln(self) -> Self;
72    /// Generic exponential for [`Float`]
73    fn exp(self) -> Self;
74
75    /// Use a primitive cast to convert a usize to the [`Float`]
76    fn usize_as_self(a: usize) -> Self;
77}
78
79macro_rules! impl_float {
80    {$($ty:ty),* } => {
81        $(
82        impl Float for $ty {
83            const MIN: Self = <$ty>::MIN;
84            const MIN_POSITIVE: Self = <$ty>::MIN_POSITIVE;
85            const MAX: Self = <$ty>::MAX;
86            const ZERO: Self = 0.0;
87            const ONE: Self = 1.0;
88            const INFINITY: Self = <$ty>::INFINITY;
89
90            #[inline]
91            fn abs(self) -> Self {
92                self.abs()
93            }
94
95            #[inline]
96            fn min(self, other: Self) -> Self {
97                self.min(other)
98            }
99
100            #[inline]
101            fn log2(self) -> Self {
102                self.log2()
103            }
104
105            #[inline]
106            fn is_nan(self) -> bool {
107                self.is_nan()
108            }
109
110            #[inline]
111            fn is_infinite(self) -> bool {
112                self.is_infinite()
113            }
114
115            #[inline]
116            fn ln(self) -> Self {
117                self.ln()
118            }
119
120            #[inline]
121            fn exp(self) -> Self {
122                self.exp()
123            }
124
125            #[inline]
126            #[allow(clippy::cast_precision_loss)]
127            fn usize_as_self(a: usize) -> $ty {
128                a as $ty
129            }
130        } )*
131
132     }
133}
134
135impl_float!(f32, f64);
136
137pub(crate) trait NearlyEqual<T> {
138    fn nearly_equal(self, b: Self, eps: T) -> bool;
139}
140
141impl<T> NearlyEqual<T> for T
142where
143    T: Float,
144{
145    #[inline]
146    /// Tests if two floating points are approximately equal within an epsilon.
147    /// Port courtesy of <https://floating-point-gui.de/errors/comparison/>
148    fn nearly_equal(self, b: Self, eps: T) -> bool {
149        let a = self;
150        let abs_a = a.abs();
151        let abs_b = b.abs();
152        let diff = (a - b).abs();
153        let zero = Self::ZERO;
154
155        if a == b {
156            // shortcut, handles infinities
157            true
158        } else if a == zero || b == zero || (abs_a + abs_b < Self::MIN_POSITIVE) {
159            // a or b is zero or both are extremely close to it
160            // relative error is less meaningful here
161            diff < eps * Self::MIN_POSITIVE
162        } else {
163            // use relative error
164            diff / (abs_a + abs_b).min(Self::MAX) < eps
165        }
166    }
167}
168
169impl<T> NearlyEqual<T> for Option<T>
170where
171    T: Float,
172{
173    #[inline]
174    fn nearly_equal(self, b: Self, eps: T) -> bool {
175        match (self, b) {
176            (Some(x), Some(y)) => x.nearly_equal(y, eps),
177            (None, None) => true,
178            _ => false,
179        }
180    }
181}
182
183/// Utlity trait for mapping floats to [`Option`].
184pub(crate) trait MapFloat: Float {
185    #[inline]
186    fn into_option(self) -> Option<Self> {
187        if self.is_infinite() || self.is_nan() {
188            None
189        } else {
190            Some(self)
191        }
192    }
193}
194
195impl MapFloat for f64 {}
196impl MapFloat for f32 {}
197
198#[macro_export]
199macro_rules! assert_fp_eq {
200    ($a:expr, $b:expr) => {
201        assert_fp_eq!($a, $b, 1e-8); // Default epsilon
202    };
203    ($a:expr, $expected:expr, $epsilon:expr) => {
204        let found = $a;
205        assert!(
206            $crate::math::NearlyEqual::nearly_equal(found, $expected, $epsilon),
207            "assertion failed: `(found ≈ expected)`\n left:\t`{:?}`,\n right:\t`{:?}`,\n eps:\t`{}`",
208            found,
209            $expected,
210            $epsilon
211        );
212    };
213}
zoe/math/float.rs

zoe/math/
float.rs
