MathUtils 1.0.0

A high-fidelity Rust mirror of the JavaScript (ECMAScript) Math built-in object.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117

//! IEEE 754-2008 binary16 (half-precision float) conversion for `Math.f16round`.

/// Rounds an `f64` value to the nearest 16-bit half-precision floating-point number (binary16),
/// returning the result as an `f64`.
///
/// This implementation corresponds to ECMAScript 2024 `Math.f16round(x)`.
///
/// # Examples
///
/// ```
/// use MathUtils::Math;
///
/// assert_eq!(Math::f16round(1.337), 1.3369140625);
/// assert_eq!(Math::f16round(std::f64::consts::PI), 3.140625);
/// assert_eq!(Math::f16round(65504.0), 65504.0);
/// assert_eq!(Math::f16round(65520.0), f64::INFINITY);
/// ```
#[inline]
pub fn f16round(x: f64) -> f64 {
    if x.is_nan() {
        return f64::NAN;
    }
    if x == 0.0 {
        return x; // Preserves +0.0 and -0.0
    }
    if !x.is_finite() {
        return x; // Preserves +Infinity and -Infinity
    }

    let bits = x.to_bits();
    let sign = (bits >> 63) & 1;
    let exp = ((bits >> 52) & 0x7ff) as i64;
    let frac = bits & 0x000f_ffff_ffff_ffff;

    // f64 exponent bias is 1023, f16 exponent bias is 15
    let new_exp = exp - 1023 + 15;

    let sign_multiplier = if sign == 1 { -1.0 } else { 1.0 };

    if new_exp >= 31 {
        // Exponent overflows f16 maximum normal range -> Infinity
        return if sign == 1 {
            f64::NEG_INFINITY
        } else {
            f64::INFINITY
        };
    }

    if new_exp <= 0 {
        // Subnormal or underflow to zero
        if new_exp < -10 {
            return if sign == 1 { -0.0 } else { 0.0 };
        }
        let full_frac = frac | 0x0010_0000_0000_0000;
        let shift = 42 + (1 - new_exp) as u64; // 52 - 10 + (1 - new_exp)
        let round_bit = (full_frac >> (shift - 1)) & 1;
        let sticky = if (full_frac & ((1 << (shift - 1)) - 1)) != 0 {
            1
        } else {
            0
        };
        let mut rounded = full_frac >> shift;
        if round_bit == 1 && (sticky == 1 || (rounded & 1) == 1) {
            rounded += 1;
        }
        let val = (rounded as f64) * 2.0f64.powi(-24);
        return sign_multiplier * val;
    }

    // Normal range
    let shift = 42u64; // 52 - 10
    let round_bit = (frac >> 41) & 1;
    let sticky = if (frac & ((1 << 41) - 1)) != 0 { 1 } else { 0 };
    let mut rounded_frac = frac >> shift;
    let mut final_exp = new_exp;
    if round_bit == 1 && (sticky == 1 || (rounded_frac & 1) == 1) {
        rounded_frac += 1;
        if rounded_frac == 1024 {
            rounded_frac = 0;
            final_exp += 1;
            if final_exp >= 31 {
                return if sign == 1 {
                    f64::NEG_INFINITY
                } else {
                    f64::INFINITY
                };
            }
        }
    }
    let val = (1.0 + (rounded_frac as f64) / 1024.0) * 2.0f64.powi((final_exp - 15) as i32);
    sign_multiplier * val
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn test_f16round_known_values() {
        assert_eq!(f16round(0.0), 0.0);
        assert!(f16round(0.0).is_sign_positive());
        assert!(f16round(-0.0).is_sign_negative());
        assert_eq!(f16round(1.0), 1.0);
        assert_eq!(f16round(-1.0), -1.0);
        assert_eq!(f16round(1.337), 1.3369140625);
        assert_eq!(f16round(std::f64::consts::PI), 3.140625);
        assert_eq!(f16round(65504.0), 65504.0);
        assert_eq!(f16round(65505.0), 65504.0);
        assert_eq!(f16round(65519.0), 65504.0);
        assert_eq!(f16round(65520.0), f64::INFINITY);
        assert_eq!(f16round(-65520.0), f64::NEG_INFINITY);
        assert_eq!(f16round(1e-8), 0.0);
        assert!(f16round(f64::NAN).is_nan());
        assert_eq!(f16round(f64::INFINITY), f64::INFINITY);
        assert_eq!(f16round(f64::NEG_INFINITY), f64::NEG_INFINITY);
    }
}