rustpython_common/
float_ops.rs

1use malachite_bigint::{BigInt, ToBigInt};
2use num_traits::{Float, Signed, ToPrimitive, Zero};
3use std::f64;
4
5pub fn ufrexp(value: f64) -> (f64, i32) {
6    if 0.0 == value {
7        (0.0, 0i32)
8    } else {
9        let bits = value.to_bits();
10        let exponent: i32 = ((bits >> 52) & 0x7ff) as i32 - 1022;
11        let mantissa_bits = bits & (0x000f_ffff_ffff_ffff) | (1022 << 52);
12        (f64::from_bits(mantissa_bits), exponent)
13    }
14}
15
16/// Equate an integer to a float.
17///
18/// Returns true if and only if, when converted to each others types, both are equal.
19///
20/// # Examples
21///
22/// ```
23/// use malachite_bigint::BigInt;
24/// use rustpython_common::float_ops::eq_int;
25/// let a = 1.0f64;
26/// let b = BigInt::from(1);
27/// let c = 2.0f64;
28/// assert!(eq_int(a, &b));
29/// assert!(!eq_int(c, &b));
30/// ```
31///
32pub fn eq_int(value: f64, other: &BigInt) -> bool {
33    if let (Some(self_int), Some(other_float)) = (value.to_bigint(), other.to_f64()) {
34        value == other_float && self_int == *other
35    } else {
36        false
37    }
38}
39
40pub fn lt_int(value: f64, other_int: &BigInt) -> bool {
41    match (value.to_bigint(), other_int.to_f64()) {
42        (Some(self_int), Some(other_float)) => value < other_float || self_int < *other_int,
43        // finite float, other_int too big for float,
44        // the result depends only on other_int’s sign
45        (Some(_), None) => other_int.is_positive(),
46        // infinite float must be bigger or lower than any int, depending on its sign
47        _ if value.is_infinite() => value.is_sign_negative(),
48        // NaN, always false
49        _ => false,
50    }
51}
52
53pub fn gt_int(value: f64, other_int: &BigInt) -> bool {
54    match (value.to_bigint(), other_int.to_f64()) {
55        (Some(self_int), Some(other_float)) => value > other_float || self_int > *other_int,
56        // finite float, other_int too big for float,
57        // the result depends only on other_int’s sign
58        (Some(_), None) => other_int.is_negative(),
59        // infinite float must be bigger or lower than any int, depending on its sign
60        _ if value.is_infinite() => value.is_sign_positive(),
61        // NaN, always false
62        _ => false,
63    }
64}
65
66pub fn div(v1: f64, v2: f64) -> Option<f64> {
67    if v2 != 0.0 {
68        Some(v1 / v2)
69    } else {
70        None
71    }
72}
73
74pub fn mod_(v1: f64, v2: f64) -> Option<f64> {
75    if v2 != 0.0 {
76        let val = v1 % v2;
77        match (v1.signum() as i32, v2.signum() as i32) {
78            (1, 1) | (-1, -1) => Some(val),
79            _ if (v1 == 0.0) || (v1.abs() == v2.abs()) => Some(val.copysign(v2)),
80            _ => Some((val + v2).copysign(v2)),
81        }
82    } else {
83        None
84    }
85}
86
87pub fn floordiv(v1: f64, v2: f64) -> Option<f64> {
88    if v2 != 0.0 {
89        Some((v1 / v2).floor())
90    } else {
91        None
92    }
93}
94
95pub fn divmod(v1: f64, v2: f64) -> Option<(f64, f64)> {
96    if v2 != 0.0 {
97        let mut m = v1 % v2;
98        let mut d = (v1 - m) / v2;
99        if v2.is_sign_negative() != m.is_sign_negative() {
100            m += v2;
101            d -= 1.0;
102        }
103        Some((d, m))
104    } else {
105        None
106    }
107}
108
109// nextafter algorithm based off of https://gitlab.com/bronsonbdevost/next_afterf
110#[allow(clippy::float_cmp)]
111pub fn nextafter(x: f64, y: f64) -> f64 {
112    if x == y {
113        y
114    } else if x.is_nan() || y.is_nan() {
115        f64::NAN
116    } else if x >= f64::INFINITY {
117        f64::MAX
118    } else if x <= f64::NEG_INFINITY {
119        f64::MIN
120    } else if x == 0.0 {
121        f64::from_bits(1).copysign(y)
122    } else {
123        // next x after 0 if y is farther from 0 than x, otherwise next towards 0
124        // the sign is a separate bit in floats, so bits+1 moves away from 0 no matter the float
125        let b = x.to_bits();
126        let bits = if (y > x) == (x > 0.0) { b + 1 } else { b - 1 };
127        let ret = f64::from_bits(bits);
128        if ret == 0.0 {
129            ret.copysign(x)
130        } else {
131            ret
132        }
133    }
134}
135
136pub fn ulp(x: f64) -> f64 {
137    if x.is_nan() {
138        return x;
139    }
140    let x = x.abs();
141    let x2 = nextafter(x, f64::INFINITY);
142    if x2.is_infinite() {
143        // special case: x is the largest positive representable float
144        let x2 = nextafter(x, f64::NEG_INFINITY);
145        x - x2
146    } else {
147        x2 - x
148    }
149}
150
151pub fn round_float_digits(x: f64, ndigits: i32) -> Option<f64> {
152    let float = if ndigits.is_zero() {
153        let fract = x.fract();
154        if (fract.abs() - 0.5).abs() < f64::EPSILON {
155            if x.trunc() % 2.0 == 0.0 {
156                x - fract
157            } else {
158                x + fract
159            }
160        } else {
161            x.round()
162        }
163    } else {
164        const NDIGITS_MAX: i32 =
165            ((f64::MANTISSA_DIGITS as i32 - f64::MIN_EXP) as f64 * f64::consts::LOG10_2) as i32;
166        const NDIGITS_MIN: i32 = -(((f64::MAX_EXP + 1) as f64 * f64::consts::LOG10_2) as i32);
167        if ndigits > NDIGITS_MAX {
168            x
169        } else if ndigits < NDIGITS_MIN {
170            0.0f64.copysign(x)
171        } else {
172            let (y, pow1, pow2) = if ndigits >= 0 {
173                // according to cpython: pow1 and pow2 are each safe from overflow, but
174                //                       pow1*pow2 ~= pow(10.0, ndigits) might overflow
175                let (pow1, pow2) = if ndigits > 22 {
176                    (10.0.powf((ndigits - 22) as f64), 1e22)
177                } else {
178                    (10.0.powf(ndigits as f64), 1.0)
179                };
180                let y = (x * pow1) * pow2;
181                if !y.is_finite() {
182                    return Some(x);
183                }
184                (y, pow1, Some(pow2))
185            } else {
186                let pow1 = 10.0.powf((-ndigits) as f64);
187                (x / pow1, pow1, None)
188            };
189            let z = y.round();
190            #[allow(clippy::float_cmp)]
191            let z = if (y - z).abs() == 0.5 {
192                2.0 * (y / 2.0).round()
193            } else {
194                z
195            };
196            let z = if let Some(pow2) = pow2 {
197                // ndigits >= 0
198                (z / pow2) / pow1
199            } else {
200                z * pow1
201            };
202
203            if !z.is_finite() {
204                // overflow
205                return None;
206            }
207
208            z
209        }
210    };
211    Some(float)
212}
rustpython_common/float_ops.rs

rustpython_common/
float_ops.rs
