1 #[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
remquo(mut x: f64, mut y: f64) -> (f64, i32)2 pub fn remquo(mut x: f64, mut y: f64) -> (f64, i32) {
3     let ux: u64 = x.to_bits();
4     let mut uy: u64 = y.to_bits();
5     let mut ex = ((ux >> 52) & 0x7ff) as i32;
6     let mut ey = ((uy >> 52) & 0x7ff) as i32;
7     let sx = (ux >> 63) != 0;
8     let sy = (uy >> 63) != 0;
9     let mut q: u32;
10     let mut i: u64;
11     let mut uxi: u64 = ux;
12 
13     if (uy << 1) == 0 || y.is_nan() || ex == 0x7ff {
14         return ((x * y) / (x * y), 0);
15     }
16     if (ux << 1) == 0 {
17         return (x, 0);
18     }
19 
20     /* normalize x and y */
21     if ex == 0 {
22         i = uxi << 12;
23         while (i >> 63) == 0 {
24             ex -= 1;
25             i <<= 1;
26         }
27         uxi <<= -ex + 1;
28     } else {
29         uxi &= (!0) >> 12;
30         uxi |= 1 << 52;
31     }
32     if ey == 0 {
33         i = uy << 12;
34         while (i >> 63) == 0 {
35             ey -= 1;
36             i <<= 1;
37         }
38         uy <<= -ey + 1;
39     } else {
40         uy &= (!0) >> 12;
41         uy |= 1 << 52;
42     }
43 
44     q = 0;
45 
46     if ex + 1 != ey {
47         if ex < ey {
48             return (x, 0);
49         }
50         /* x mod y */
51         while ex > ey {
52             i = uxi.wrapping_sub(uy);
53             if (i >> 63) == 0 {
54                 uxi = i;
55                 q += 1;
56             }
57             uxi <<= 1;
58             q <<= 1;
59             ex -= 1;
60         }
61         i = uxi.wrapping_sub(uy);
62         if (i >> 63) == 0 {
63             uxi = i;
64             q += 1;
65         }
66         if uxi == 0 {
67             ex = -60;
68         } else {
69             while (uxi >> 52) == 0 {
70                 uxi <<= 1;
71                 ex -= 1;
72             }
73         }
74     }
75 
76     /* scale result and decide between |x| and |x|-|y| */
77     if ex > 0 {
78         uxi -= 1 << 52;
79         uxi |= (ex as u64) << 52;
80     } else {
81         uxi >>= -ex + 1;
82     }
83     x = f64::from_bits(uxi);
84     if sy {
85         y = -y;
86     }
87     if ex == ey || (ex + 1 == ey && (2.0 * x > y || (2.0 * x == y && (q % 2) != 0))) {
88         x -= y;
89         // TODO: this matches musl behavior, but it is incorrect
90         q = q.wrapping_add(1);
91     }
92     q &= 0x7fffffff;
93     let quo = if sx ^ sy { -(q as i32) } else { q as i32 };
94     if sx {
95         (-x, quo)
96     } else {
97         (x, quo)
98     }
99 }
100 
101 #[cfg(test)]
102 mod tests {
103     use super::remquo;
104 
105     #[test]
test_q_overflow()106     fn test_q_overflow() {
107         // 0xc000000000000001, 0x04c0000000000004
108         let _ = remquo(-2.0000000000000004, 8.406091369059082e-286);
109     }
110 }
111