qiurui
/
glibc
mirror of git://sourceware.org/git/glibc.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

math: Simplify and optimize modf implementation 

Refactor the generic implementation to use math_config.h definitions,
and add an alternative one if the ABI supports truncf instructions
(gated through math-use-builtins-trunc.h).

The generic implementation generates similar code on x86_64, while
the optimization one for aarch64 (where truncf is supported as a
builtin by through frintz), the improvements are:

reciprocal-throughput           master    patch    difference
workload-0_1                    3.0595   3.0698        -0.34%
workload-1_maxint               5.1747   3.0542        40.98%
workload-maxint_maxfloat        3.4391   3.0349        11.75%
workload-integral               3.2732   3.0293         7.45%

latency                         master    patch    difference
workload-0_1                    3.5267   4.7107       -33.57%
workload-1_maxint               6.9074   4.7282        31.55%
workload-maxint_maxfloat        3.7210   4.7506       -27.67%
workload-integral               3.8634   4.8137       -24.60%

Checked on aarch64-linux-gnu and x86_64-linux-gnu.
Reviewed-by: Wilco Dijkstra  <Wilco.Dijkstra@arm.com>
This commit is contained in:
Adhemerval Zanella 2025-06-16 10:17:34 -03:00
parent 61cc9922f3
commit f165e244e4
2 changed files with 70 additions and 52 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

13
sysdeps/ieee754/dbl-64/math_config.h
View File

@ -109,6 +109,7 @@ issignaling_inline (double x)
#define BIT_WIDTH 64
#define MANTISSA_WIDTH 52
#define EXPONENT_WIDTH 11
#define EXPONENT_BIAS 1023
#define MANTISSA_MASK UINT64_C(0x000fffffffffffff)
#define EXPONENT_MASK UINT64_C(0x7ff0000000000000)
#define EXP_MANT_MASK UINT64_C(0x7fffffffffffffff)
@ -121,12 +122,24 @@ is_nan (uint64_t x)
return (x & EXP_MANT_MASK) > EXPONENT_MASK;
}
static inline bool
is_inf (uint64_t x)
{
return (x << 1) == (EXPONENT_MASK << 1);
}
static inline uint64_t
get_mantissa (uint64_t x)
{
return x & MANTISSA_MASK;
}
static inline int
get_exponent (uint64_t x)
{
return (int)((x >> MANTISSA_WIDTH & 0x7ff) - EXPONENT_BIAS);
}
/* Convert integer number X, unbiased exponent EP, and sign S to double:
result = X * 2^(EP+1 - exponent_bias)

109
sysdeps/ieee754/dbl-64/s_modf.c
View File

@ -1,63 +1,68 @@
/*
* ====================================================
* 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.
* ====================================================
*/
/* Extract signed integral and fractional values.
Copyright (C) 1993-2025 Free Software Foundation, Inc.
This file is part of the GNU C Library.
/*
* modf(double x, double *iptr)
* return fraction part of x, and return x's integral part in *iptr.
* Method:
* Bit twiddling.
*
* Exception:
* No exception.
*/
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
<https://www.gnu.org/licenses/>. */
#include <math.h>
#include <math_private.h>
#include <libm-alias-double.h>
#include <stdint.h>
static const double one = 1.0;
#include "math_config.h"
#include <math-use-builtins-trunc.h>
double
__modf(double x, double *iptr)
__modf (double x, double *iptr)
{
int64_t i0;
int32_t j0;
EXTRACT_WORDS64(i0,x);
j0 = ((i0>>52)&0x7ff)-0x3ff; /* exponent of x */
if(j0<52) { /* integer part in x */
if(j0<0) { /* |x|<1 */
/* *iptr = +-0 */
INSERT_WORDS64(*iptr,i0&UINT64_C(0x8000000000000000));
return x;
} else {
uint64_t i = UINT64_C(0x000fffffffffffff)>>j0;
if((i0&i)==0) { /* x is integral */
*iptr = x;
/* return +-0 */
INSERT_WORDS64(x,i0&UINT64_C(0x8000000000000000));
return x;
} else {
INSERT_WORDS64(*iptr,i0&(~i));
return x - *iptr;
}
}
} else { /* no fraction part */
*iptr = x*one;
/* We must handle NaNs separately. */
if (j0 == 0x400 && (i0 & UINT64_C(0xfffffffffffff)))
return x*one;
INSERT_WORDS64(x,i0&UINT64_C(0x8000000000000000)); /* return +-0 */
return x;
uint64_t t = asuint64 (x);
#if USE_TRUNC_BUILTIN
if (is_inf (t))
{
*iptr = x;
return copysign (0.0, x);
}
*iptr = trunc (x);
return copysign (x - *iptr, x);
#else
int e = get_exponent (t);
/* No fraction part. */
if (e < MANTISSA_WIDTH)
{
if (e < 0)
{
/* |x|<1 -> *iptr = +-0 */
*iptr = asdouble (t & SIGN_MASK);
return x;
}
uint64_t i = MANTISSA_MASK >> e;
if ((t & i) == 0)
{
/* x in integral, return +-0 */
*iptr = x;
return asdouble (t & SIGN_MASK);
}
*iptr = asdouble (t & ~i);
return x - *iptr;
}
/* Set invalid operation for sNaN. */
*iptr = x * 1.0;
if ((e == 0x400) && (t & MANTISSA_MASK))
return *iptr;
return asdouble (t & SIGN_MASK);
#endif
}
#ifndef __modf
libm_alias_double (__modf, modf)