Skip to content
This repository was archived by the owner on Apr 28, 2025. It is now read-only.

Commit d48187a

Browse files
tgross35tgross35
committed
Add a generic version of fmod
This can replace `fmod` and `fmodf`.
1 parent 42a8dfb commit d48187a

File tree

6 files changed

+107
-162
lines changed

6 files changed

+107
-162
lines changed

‎etc/function-definitions.json

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -437,13 +437,15 @@
437437
"fmod": {
438438
"sources": [
439439
"src/libm_helper.rs",
440-
"src/math/fmod.rs"
440+
"src/math/fmod.rs",
441+
"src/math/generic/fmod.rs"
441442
],
442443
"type": "f64"
443444
},
444445
"fmodf": {
445446
"sources": [
446-
"src/math/fmodf.rs"
447+
"src/math/fmodf.rs",
448+
"src/math/generic/fmod.rs"
447449
],
448450
"type": "f32"
449451
},

‎src/math/fmod.rs

Lines changed: 2 additions & 75 deletions
Original file line numberDiff line numberDiff line change
@@ -1,78 +1,5 @@
1+
/// Calculate the remainder of `x / y`, the precise result of `x - trunc(x / y) * y`.
12
#[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
23
pub fn fmod(x: f64, y: f64) -> f64 {
3-
let mut uxi = x.to_bits();
4-
let mut uyi = y.to_bits();
5-
let mut ex = ((uxi >> 52) & 0x7ff) as i64;
6-
let mut ey = ((uyi >> 52) & 0x7ff) as i64;
7-
let sx = uxi >> 63;
8-
let mut i;
9-
10-
if uyi << 1 == 0 || y.is_nan() || ex == 0x7ff {
11-
return (x * y) / (x * y);
12-
}
13-
if uxi << 1 <= uyi << 1 {
14-
if uxi << 1 == uyi << 1 {
15-
return 0.0 * x;
16-
}
17-
return x;
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 &= u64::MAX >> 12;
30-
uxi |= 1 << 52;
31-
}
32-
if ey == 0 {
33-
i = uyi << 12;
34-
while i >> 63 == 0 {
35-
ey -= 1;
36-
i <<= 1;
37-
}
38-
uyi <<= -ey + 1;
39-
} else {
40-
uyi &= u64::MAX >> 12;
41-
uyi |= 1 << 52;
42-
}
43-
44-
/* x mod y */
45-
while ex > ey {
46-
i = uxi.wrapping_sub(uyi);
47-
if i >> 63 == 0 {
48-
if i == 0 {
49-
return 0.0 * x;
50-
}
51-
uxi = i;
52-
}
53-
uxi <<= 1;
54-
ex -= 1;
55-
}
56-
i = uxi.wrapping_sub(uyi);
57-
if i >> 63 == 0 {
58-
if i == 0 {
59-
return 0.0 * x;
60-
}
61-
uxi = i;
62-
}
63-
while uxi >> 52 == 0 {
64-
uxi <<= 1;
65-
ex -= 1;
66-
}
67-
68-
/* scale result */
69-
if ex > 0 {
70-
uxi -= 1 << 52;
71-
uxi |= (ex as u64) << 52;
72-
} else {
73-
uxi >>= -ex + 1;
74-
}
75-
uxi |= sx << 63;
76-
77-
f64::from_bits(uxi)
4+
super::generic::fmod(x, y)
785
}

‎src/math/fmodf.rs

Lines changed: 2 additions & 85 deletions
Original file line numberDiff line numberDiff line change
@@ -1,88 +1,5 @@
1-
use core::f32;
2-
1+
/// Calculate the remainder of `x / y`, the precise result of `x - trunc(x / y) * y`.
32
#[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
43
pub fn fmodf(x: f32, y: f32) -> f32 {
5-
let mut uxi = x.to_bits();
6-
let mut uyi = y.to_bits();
7-
let mut ex = ((uxi >> 23) & 0xff) as i32;
8-
let mut ey = ((uyi >> 23) & 0xff) as i32;
9-
let sx = uxi & 0x80000000;
10-
let mut i;
11-
12-
if uyi << 1 == 0 || y.is_nan() || ex == 0xff {
13-
return (x * y) / (x * y);
14-
}
15-
16-
if uxi << 1 <= uyi << 1 {
17-
if uxi << 1 == uyi << 1 {
18-
return 0.0 * x;
19-
}
20-
21-
return x;
22-
}
23-
24-
/* normalize x and y */
25-
if ex == 0 {
26-
i = uxi << 9;
27-
while i >> 31 == 0 {
28-
ex -= 1;
29-
i <<= 1;
30-
}
31-
32-
uxi <<= -ex + 1;
33-
} else {
34-
uxi &= u32::MAX >> 9;
35-
uxi |= 1 << 23;
36-
}
37-
38-
if ey == 0 {
39-
i = uyi << 9;
40-
while i >> 31 == 0 {
41-
ey -= 1;
42-
i <<= 1;
43-
}
44-
45-
uyi <<= -ey + 1;
46-
} else {
47-
uyi &= u32::MAX >> 9;
48-
uyi |= 1 << 23;
49-
}
50-
51-
/* x mod y */
52-
while ex > ey {
53-
i = uxi.wrapping_sub(uyi);
54-
if i >> 31 == 0 {
55-
if i == 0 {
56-
return 0.0 * x;
57-
}
58-
uxi = i;
59-
}
60-
uxi <<= 1;
61-
62-
ex -= 1;
63-
}
64-
65-
i = uxi.wrapping_sub(uyi);
66-
if i >> 31 == 0 {
67-
if i == 0 {
68-
return 0.0 * x;
69-
}
70-
uxi = i;
71-
}
72-
73-
while uxi >> 23 == 0 {
74-
uxi <<= 1;
75-
ex -= 1;
76-
}
77-
78-
/* scale result up */
79-
if ex > 0 {
80-
uxi -= 1 << 23;
81-
uxi |= (ex as u32) << 23;
82-
} else {
83-
uxi >>= -ex + 1;
84-
}
85-
uxi |= sx;
86-
87-
f32::from_bits(uxi)
4+
super::generic::fmod(x, y)
885
}

‎src/math/generic/fmod.rs

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

‎src/math/generic/mod.rs

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -5,6 +5,7 @@ mod fdim;
55
mod floor;
66
mod fmax;
77
mod fmin;
8+
mod fmod;
89
mod rint;
910
mod round;
1011
mod scalbn;
@@ -18,6 +19,7 @@ pub use fdim::fdim;
1819
pub use floor::floor;
1920
pub use fmax::fmax;
2021
pub use fmin::fmin;
22+
pub use fmod::fmod;
2123
pub use rint::rint;
2224
pub use round::round;
2325
pub use scalbn::scalbn;

‎src/math/support/int_traits.rs

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -46,6 +46,7 @@ pub trait Int:
4646
+ ops::BitXorAssign
4747
+ ops::ShlAssign<i32>
4848
+ ops::ShrAssign<u32>
49+
+ ops::ShrAssign<i32>
4950
+ ops::Add<Output = Self>
5051
+ ops::Sub<Output = Self>
5152
+ ops::Mul<Output = Self>

0 commit comments

Comments
 (0)