Skip to main content

c/
math_bitmanip.rs

1// Copyright (c) 2026 Meng Zhuo
2// SPDX-License-Identifier: BSD-3-Clause
3
4#[inline]
5fn p2d(e: i32) -> f64 {
6    f64::from_bits(((0x3ff + e) as u64) << 52)
7}
8
9#[inline]
10fn p2f(e: i32) -> f32 {
11    f32::from_bits(((0x7f + e) as u32) << 23)
12}
13
14#[no_mangle]
15pub extern "C" fn fabs(x: f64) -> f64 {
16    f64::from_bits(x.to_bits() & ((-1i64 as u64) / 2))
17}
18
19#[no_mangle]
20pub extern "C" fn fabsf(x: f32) -> f32 {
21    f32::from_bits(x.to_bits() & 0x7fffffff)
22}
23
24#[no_mangle]
25pub extern "C" fn copysign(x: f64, y: f64) -> f64 {
26    let mut ux = x.to_bits();
27    ux &= (-1i64 as u64) / 2;
28    ux |= y.to_bits() & (1u64 << 63);
29    f64::from_bits(ux)
30}
31
32#[no_mangle]
33pub extern "C" fn copysignf(x: f32, y: f32) -> f32 {
34    let mut ux = x.to_bits();
35    ux &= 0x7fffffff;
36    ux |= y.to_bits() & 0x80000000;
37    f32::from_bits(ux)
38}
39
40#[no_mangle]
41pub extern "C" fn trunc(x: f64) -> f64 {
42    let ix = x.to_bits();
43    let mut e = ((ix >> 52) & 0x7ff) as i32 - 0x3ff + 12;
44
45    if e >= 52 + 12 {
46        return x;
47    }
48    if e < 12 {
49        e = 1;
50    }
51    let m = (-1i64 as u64) >> e;
52    if (ix & m) == 0 {
53        return x;
54    }
55    force_eval(x + p2d(120));
56    f64::from_bits(ix & !m)
57}
58
59#[no_mangle]
60pub extern "C" fn truncf(x: f32) -> f32 {
61    let ix = x.to_bits();
62    let mut e = ((ix >> 23) & 0xff) as i32 - 0x7f + 9;
63
64    if e >= 23 + 9 {
65        return x;
66    }
67    if e < 9 {
68        e = 1;
69    }
70    let m = (-1i32 as u32) >> e;
71    if (ix & m) == 0 {
72        return x;
73    }
74    force_eval(x + p2f(120));
75    f32::from_bits(ix & !m)
76}
77
78#[no_mangle]
79pub extern "C" fn floor(x: f64) -> f64 {
80    let ix = x.to_bits();
81    let e = ((ix >> 52) & 0x7ff) as i32;
82
83    if e >= 0x3ff + 52 || x == 0.0 {
84        return x;
85    }
86    let y = if ix >> 63 != 0 {
87        x - p2d(52) + p2d(52) - x
88    } else {
89        x + p2d(52) - p2d(52) - x
90    };
91    if e <= 0x3ff - 1 {
92        force_eval(y);
93        return if ix >> 63 != 0 { -1.0 } else { 0.0 };
94    }
95    if y > 0.0 {
96        x + y - 1.0
97    } else {
98        x + y
99    }
100}
101
102#[no_mangle]
103pub extern "C" fn floorf(x: f32) -> f32 {
104    let mut ix = x.to_bits();
105    let e = ((ix >> 23) & 0xff) as i32 - 0x7f;
106
107    if e >= 23 {
108        return x;
109    }
110    if e >= 0 {
111        let m = 0x007fffff >> e;
112        if (ix & m) == 0 {
113            return x;
114        }
115        force_eval(x + p2f(120));
116        if ix >> 31 != 0 {
117            ix += m;
118        }
119        ix &= !m;
120    } else {
121        force_eval(x + p2f(120));
122        if ix >> 31 == 0 {
123            ix = 0;
124        } else if ix << 1 != 0 {
125            return -1.0f32;
126        }
127    }
128    f32::from_bits(ix)
129}
130
131#[no_mangle]
132pub extern "C" fn ceil(x: f64) -> f64 {
133    let ix = x.to_bits();
134    let e = ((ix >> 52) & 0x7ff) as i32;
135
136    if e >= 0x3ff + 52 || x == 0.0 {
137        return x;
138    }
139    let y = if ix >> 63 != 0 {
140        x - p2d(52) + p2d(52) - x
141    } else {
142        x + p2d(52) - p2d(52) - x
143    };
144    if e <= 0x3ff - 1 {
145        force_eval(y);
146        return if ix >> 63 != 0 { -0.0 } else { 1.0 };
147    }
148    if y < 0.0 {
149        x + y + 1.0
150    } else {
151        x + y
152    }
153}
154
155#[no_mangle]
156pub extern "C" fn ceilf(x: f32) -> f32 {
157    let mut ix = x.to_bits();
158    let e = ((ix >> 23) & 0xff) as i32 - 0x7f;
159
160    if e >= 23 {
161        return x;
162    }
163    if e >= 0 {
164        let m = 0x007fffff >> e;
165        if (ix & m) == 0 {
166            return x;
167        }
168        force_eval(x + p2f(120));
169        if ix >> 31 == 0 {
170            ix += m;
171        }
172        ix &= !m;
173    } else {
174        force_eval(x + p2f(120));
175        if ix >> 31 != 0 {
176            return -0.0f32;
177        } else if ix << 1 != 0 {
178            return 1.0f32;
179        }
180    }
181    f32::from_bits(ix)
182}
183
184#[no_mangle]
185pub extern "C" fn round(x: f64) -> f64 {
186    let ix = x.to_bits();
187    let e = ((ix >> 52) & 0x7ff) as i32;
188    let mut ax = x;
189
190    if e >= 0x3ff + 52 {
191        return x;
192    }
193    let neg = ix >> 63 != 0;
194    if neg {
195        ax = -ax;
196    }
197    if e < 0x3ff - 1 {
198        force_eval(ax + p2d(52));
199        return 0.0 * x;
200    }
201    let y = ax + p2d(52) - p2d(52) - ax;
202    let y = if y > 0.5 {
203        y + ax - 1.0
204    } else if y <= -0.5 {
205        y + ax + 1.0
206    } else {
207        y + ax
208    };
209    if neg { -y } else { y }
210}
211
212#[no_mangle]
213pub extern "C" fn roundf(x: f32) -> f32 {
214    let ix = x.to_bits();
215    let e = ((ix >> 23) & 0xff) as i32;
216    let mut ax = x;
217
218    if e >= 0x7f + 23 {
219        return x;
220    }
221    let neg = ix >> 31 != 0;
222    if neg {
223        ax = -ax;
224    }
225    if e < 0x7f - 1 {
226        force_eval(ax + p2f(23));
227        return 0.0f32 * x;
228    }
229    let y = ax + p2f(23) - p2f(23) - ax;
230    let y = if y > 0.5f32 {
231        y + ax - 1.0f32
232    } else if y <= -0.5f32 {
233        y + ax + 1.0f32
234    } else {
235        y + ax
236    };
237    if neg { -y } else { y }
238}
239
240#[no_mangle]
241pub extern "C" fn scalbn(x: f64, n: c_int) -> f64 {
242    let mut y = x;
243    let mut n = n;
244
245    if n > 1023 {
246        y *= p2d(1023);
247        n -= 1023;
248        if n > 1023 {
249            y *= p2d(1023);
250            n -= 1023;
251            if n > 1023 {
252                n = 1023;
253            }
254        }
255    } else if n < -1022 {
256        y *= p2d(-1022) * p2d(53);
257        n += 1022 - 53;
258        if n < -1022 {
259            y *= p2d(-1022) * p2d(53);
260            n += 1022 - 53;
261            if n < -1022 {
262                n = -1022;
263            }
264        }
265    }
266    f64::from_bits(((0x3ff + n) as u64) << 52) * y
267}
268
269#[no_mangle]
270pub extern "C" fn scalbnf(x: f32, n: c_int) -> f32 {
271    let mut y = x;
272    let mut n = n;
273
274    if n > 127 {
275        y *= p2f(127);
276        n -= 127;
277        if n > 127 {
278            y *= p2f(127);
279            n -= 127;
280            if n > 127 {
281                n = 127;
282            }
283        }
284    } else if n < -126 {
285        y *= p2f(-126) * p2f(24);
286        n += 126 - 24;
287        if n < -126 {
288            y *= p2f(-126) * p2f(24);
289            n += 126 - 24;
290            if n < -126 {
291                n = -126;
292            }
293        }
294    }
295    f32::from_bits(((0x7f + n) as u32) << 23) * y
296}
297
298#[no_mangle]
299pub extern "C" fn scalbln(x: f64, n: c_long) -> f64 {
300    let n = if n > c_int::MAX as c_long {
301        c_int::MAX
302    } else if n < c_int::MIN as c_long {
303        c_int::MIN
304    } else {
305        n as c_int
306    };
307    scalbn(x, n)
308}
309
310#[no_mangle]
311pub extern "C" fn scalblnf(x: f32, n: c_long) -> f32 {
312    let n = if n > c_int::MAX as c_long {
313        c_int::MAX
314    } else if n < c_int::MIN as c_long {
315        c_int::MIN
316    } else {
317        n as c_int
318    };
319    scalbnf(x, n)
320}
321
322#[no_mangle]
323pub extern "C" fn ldexp(x: f64, n: c_int) -> f64 {
324    scalbn(x, n)
325}
326
327#[no_mangle]
328pub extern "C" fn ldexpf(x: f32, n: c_int) -> f32 {
329    scalbnf(x, n)
330}
331
332#[no_mangle]
333pub unsafe extern "C" fn frexp(x: f64, e: *mut c_int) -> f64 {
334    let mut ix = x.to_bits();
335    let ee = ((ix >> 52) & 0x7ff) as i32;
336
337    if ee == 0 {
338        if x != 0.0 {
339            let r = frexp(x * p2d(64), e);
340            if !e.is_null() {
341                *e -= 64;
342            }
343            return r;
344        }
345        if !e.is_null() {
346            *e = 0;
347        }
348        return x;
349    }
350    if ee == 0x7ff {
351        if !e.is_null() {
352            *e = 0;
353        }
354        return x;
355    }
356    if !e.is_null() {
357        *e = ee - 0x3fe;
358    }
359    ix &= 0x800fffffffffffffu64;
360    ix |= 0x3fe0000000000000u64;
361    f64::from_bits(ix)
362}
363
364#[no_mangle]
365pub unsafe extern "C" fn frexpf(x: f32, e: *mut c_int) -> f32 {
366    let mut ix = x.to_bits();
367    let ee = ((ix >> 23) & 0xff) as i32;
368
369    if ee == 0 {
370        if x != 0.0f32 {
371            let r = frexpf(x * p2f(64), e);
372            if !e.is_null() {
373                *e -= 64;
374            }
375            return r;
376        }
377        if !e.is_null() {
378            *e = 0;
379        }
380        return x;
381    }
382    if ee == 0xff {
383        if !e.is_null() {
384            *e = 0;
385        }
386        return x;
387    }
388    if !e.is_null() {
389        *e = ee - 0x7e;
390    }
391    ix &= 0x807fffffu32;
392    ix |= 0x3f000000u32;
393    f32::from_bits(ix)
394}
395
396#[no_mangle]
397pub unsafe extern "C" fn modf(x: f64, iptr: *mut f64) -> f64 {
398    let mut ix = x.to_bits();
399    let e = ((ix >> 52) & 0x7ff) as i32 - 0x3ff;
400
401    if e >= 52 {
402        if !iptr.is_null() {
403            *iptr = x;
404        }
405        if e == 0x400 && (ix << 12) != 0 {
406            return x;
407        }
408        return f64::from_bits(ix & (1u64 << 63));
409    }
410    if e < 0 {
411        if !iptr.is_null() {
412            *iptr = f64::from_bits(ix & (1u64 << 63));
413        }
414        return x;
415    }
416    let mask = (-1i64 as u64) >> 12 >> e;
417    if (ix & mask) == 0 {
418        if !iptr.is_null() {
419            *iptr = x;
420        }
421        return f64::from_bits(ix & (1u64 << 63));
422    }
423    ix &= !mask;
424    let i = f64::from_bits(ix);
425    if !iptr.is_null() {
426        *iptr = i;
427    }
428    x - i
429}
430
431#[no_mangle]
432pub unsafe extern "C" fn modff(x: f32, iptr: *mut f32) -> f32 {
433    let mut ix = x.to_bits();
434    let e = ((ix >> 23) & 0xff) as i32 - 0x7f;
435
436    if e >= 23 {
437        if !iptr.is_null() {
438            *iptr = x;
439        }
440        if e == 0x80 && (ix << 9) != 0 {
441            return x;
442        }
443        return f32::from_bits(ix & 0x80000000u32);
444    }
445    if e < 0 {
446        if !iptr.is_null() {
447            *iptr = f32::from_bits(ix & 0x80000000u32);
448        }
449        return x;
450    }
451    let mask = 0x007fffffu32 >> e;
452    if (ix & mask) == 0 {
453        if !iptr.is_null() {
454            *iptr = x;
455        }
456        return f32::from_bits(ix & 0x80000000u32);
457    }
458    ix &= !mask;
459    let i = f32::from_bits(ix);
460    if !iptr.is_null() {
461        *iptr = i;
462    }
463    x - i
464}