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