1/* SPDX-License-Identifier: MIT */
2/* origin: musl src/math/fmod.c. Ported to generic Rust algorithm in 2025, TG. */
3
4use super::super::{CastFrom, Float, Int, MinInt};
5
6#[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
7pub fn fmod<F: Float>(x: F, y: F) -> F {
8 let zero = F::Int::ZERO;
9 let one = F::Int::ONE;
10 let mut ix = x.to_bits();
11 let mut iy = y.to_bits();
12 let mut ex = x.ex().signed();
13 let mut ey = y.ex().signed();
14 let sx = ix & F::SIGN_MASK;
15
16 if iy << 1 == zero || y.is_nan() || ex == F::EXP_SAT as i32 {
17 return (x * y) / (x * y);
18 }
19
20 if ix << 1 <= iy << 1 {
21 if ix << 1 == iy << 1 {
22 return F::ZERO * x;
23 }
24 return x;
25 }
26
27 /* normalize x and y */
28 if ex == 0 {
29 let i = ix << F::EXP_BITS;
30 ex -= i.leading_zeros() as i32;
31 ix <<= -ex + 1;
32 } else {
33 ix &= F::Int::MAX >> F::EXP_BITS;
34 ix |= one << F::SIG_BITS;
35 }
36
37 if ey == 0 {
38 let i = iy << F::EXP_BITS;
39 ey -= i.leading_zeros() as i32;
40 iy <<= -ey + 1;
41 } else {
42 iy &= F::Int::MAX >> F::EXP_BITS;
43 iy |= one << F::SIG_BITS;
44 }
45
46 /* x mod y */
47 while ex > ey {
48 let i = ix.wrapping_sub(iy);
49 if i >> (F::BITS - 1) == zero {
50 if i == zero {
51 return F::ZERO * x;
52 }
53 ix = i;
54 }
55
56 ix <<= 1;
57 ex -= 1;
58 }
59
60 let i = ix.wrapping_sub(iy);
61 if i >> (F::BITS - 1) == zero {
62 if i == zero {
63 return F::ZERO * x;
64 }
65
66 ix = i;
67 }
68
69 let shift = ix.leading_zeros().saturating_sub(F::EXP_BITS);
70 ix <<= shift;
71 ex -= shift as i32;
72
73 /* scale result */
74 if ex > 0 {
75 ix -= one << F::SIG_BITS;
76 ix |= F::Int::cast_from(ex) << F::SIG_BITS;
77 } else {
78 ix >>= -ex + 1;
79 }
80
81 ix |= sx;
82
83 F::from_bits(ix)
84}
85