1 use core::u64;
2 
3 #[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
fmod(x: f64, y: f64) -> f644 pub 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