From 4d2582150e4995c4ff0c4e9f678a4fed02830513 Mon Sep 17 00:00:00 2001 From: Osama Abdelkader Date: Thu, 23 Oct 2025 18:06:28 +0300 Subject: [PATCH] math: Optimize frexpf (binary32) with fast path for normal numbers Add fast path optimization for frexpf using a single unsigned comparison to identify normal floating-point numbers and return immediately via arithmetic on the bit representation. The implementation uses asuint()/asfloat() from math_config.h and arithmetic operations to adjust the exponent, which generates better code than bit masking on ARM and RISC-V architectures. For subnormals, stdc_leading_zeros provides faster normalization than the traditional multiply approach. The zero/infinity/NaN check is simplified to (int32_t)(hx << 1) <= 0, which is more efficient than separate comparisons. Benchmark results on Intel Core i9-13900H (13th Gen): Baseline: 5.858 ns/op Optimized: 4.003 ns/op Speedup: 1.46x (31.7% faster) Zero: 3.580 ns/op (fast path) Denormal: 5.597 ns/op (slower, rare case) Signed-off-by: Osama Abdelkader Reviewed-by: Adhemerval Zanella --- sysdeps/ieee754/flt-32/s_frexpf.c | 82 +++++++++++++++++-------------- 1 file changed, 46 insertions(+), 36 deletions(-) diff --git a/sysdeps/ieee754/flt-32/s_frexpf.c b/sysdeps/ieee754/flt-32/s_frexpf.c index 59fef66a6b8..f424e52199f 100644 --- a/sysdeps/ieee754/flt-32/s_frexpf.c +++ b/sysdeps/ieee754/flt-32/s_frexpf.c @@ -1,44 +1,54 @@ -/* s_frexpf.c -- float version of s_frexp.c. - */ - -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunPro, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ - -#if defined(LIBM_SCCS) && !defined(lint) -static char rcsid[] = "$NetBSD: s_frexpf.c,v 1.5 1995/05/10 20:47:26 jtc Exp $"; -#endif +/* Optimized frexp implementation for binary32. + Copyright (C) 2025 Free Software Foundation, Inc. + This file is part of the GNU C Library. + + The GNU C Library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 2.1 of the License, or (at your option) any later version. + + The GNU C Library is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public + License along with the GNU C Library; if not, see + . */ #include #include +#include +#include "math_config.h" #include -static const float -two25 = 3.3554432000e+07; /* 0x4c000000 */ - -float __frexpf(float x, int *eptr) +float +__frexpf (float x, int *eptr) { - int32_t hx,ix; - GET_FLOAT_WORD(hx,x); - ix = 0x7fffffff&hx; - *eptr = 0; - if(ix>=0x7f800000||(ix==0)) return x + x; /* 0,inf,nan */ - if (ix<0x00800000) { /* subnormal */ - x *= two25; - GET_FLOAT_WORD(hx,x); - ix = hx&0x7fffffff; - *eptr = -25; - } - *eptr += (ix>>23)-126; - hx = (hx&0x807fffff)|0x3f000000; - SET_FLOAT_WORD(x,hx); - return x; + uint32_t hx = asuint (x); + uint32_t ex = (hx >> MANTISSA_WIDTH) & 0xff; + + /* Fast path for normal numbers. */ + if (__glibc_likely ((ex - 1) < 0xfe)) + { + int e = ex - EXPONENT_BIAS + 1; + *eptr = e; + return asfloat (hx - (e << MANTISSA_WIDTH)); + } + + /* Handle zero, infinity, and NaN. */ + if (__glibc_likely ((int32_t) (hx << 1) <= 0)) + { + *eptr = 0; + return x + x; + } + + /* Subnormal. */ + uint32_t sign = hx & SIGN_MASK; + int lz = stdc_leading_zeros (hx << (32 - MANTISSA_WIDTH - 1)); + hx <<= lz; + *eptr = -(EXPONENT_BIAS - 2) - lz; + return asfloat ((hx & MANTISSA_MASK) | sign + | (((uint32_t) (EXPONENT_BIAS - 1)) << MANTISSA_WIDTH)); } libm_alias_float (__frexp, frexp) -- 2.47.3