1use core::u64;
2
3#[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
4pub fn fmod(x: f64, y: f64) -> f64 {
5 let mut uxi = x.to_bits();
6 let mut uyi = y.to_bits();
7 let mut ex = (uxi >> 52 & 0x7ff) as i64;
8 let mut ey = (uyi >> 52 & 0x7ff) as i64;
9 let sx = uxi >> 63;
10 let mut i;
11
12 if uyi << 1 == 0 || y.is_nan() || ex == 0x7ff {
13 return (x * y) / (x * y);
14 }
15 if uxi << 1 <= uyi << 1 {
16 if uxi << 1 == uyi << 1 {
17 return 0.0 * x;
18 }
19 return x;
20 }
21
22 /* normalize x and y */
23 if ex == 0 {
24 i = uxi << 12;
25 while i >> 63 == 0 {
26 ex -= 1;
27 i <<= 1;
28 }
29 uxi <<= -ex + 1;
30 } else {
31 uxi &= u64::MAX >> 12;
32 uxi |= 1 << 52;
33 }
34 if ey == 0 {
35 i = uyi << 12;
36 while i >> 63 == 0 {
37 ey -= 1;
38 i <<= 1;
39 }
40 uyi <<= -ey + 1;
41 } else {
42 uyi &= u64::MAX >> 12;
43 uyi |= 1 << 52;
44 }
45
46 /* x mod y */
47 while ex > ey {
48 i = uxi.wrapping_sub(uyi);
49 if i >> 63 == 0 {
50 if i == 0 {
51 return 0.0 * x;
52 }
53 uxi = i;
54 }
55 uxi <<= 1;
56 ex -= 1;
57 }
58 i = uxi.wrapping_sub(uyi);
59 if i >> 63 == 0 {
60 if i == 0 {
61 return 0.0 * x;
62 }
63 uxi = i;
64 }
65 while uxi >> 52 == 0 {
66 uxi <<= 1;
67 ex -= 1;
68 }
69
70 /* scale result */
71 if ex > 0 {
72 uxi -= 1 << 52;
73 uxi |= (ex as u64) << 52;
74 } else {
75 uxi >>= -ex + 1;
76 }
77 uxi |= (sx as u64) << 63;
78
79 f64::from_bits(uxi)
80}
81