glibc/math/s_exp2m1_template.c

66 lines
1.9 KiB
C
Raw Normal View History

Implement C23 exp2m1, exp10m1 C23 adds various <math.h> function families originally defined in TS 18661-4. Add the exp2m1 and exp10m1 functions (exp2(x)-1 and exp10(x)-1, like expm1). As with other such functions, these use type-generic templates that could be replaced with faster and more accurate type-specific implementations in future. Test inputs are copied from those for expm1, plus some additions close to the overflow threshold (copied from exp2 and exp10) and also some near the underflow threshold. exp2m1 has the unusual property of having an input (M_MAX_EXP) where whether the function overflows (under IEEE semantics) depends on the rounding mode. Although these could reasonably be XFAILed in the testsuite (as we do in some cases for arguments very close to a function's overflow threshold when an error of a few ulps in the implementation can result in the implementation not agreeing with an ideal one on whether overflow takes place - the testsuite isn't smart enough to handle this automatically), since these functions aren't required to be correctly rounding, I made the implementation check for and handle this case specially. The Makefile ordering expected by lint-makefiles for the new functions is a bit peculiar, but I implemented it in this patch so that the test passes; I don't know why log2 also needed moving in one Makefile variable setting when it didn't in my previous patches, but the failure showed a different place was expected for that function as well. The powerpc64le IFUNC setup seems not to be as self-contained as one might hope; it shouldn't be necessary to add IFUNCs for new functions such as these simply to get them building, but without setting up IFUNCs for the new functions, there were undefined references to __GI___expm1f128 (that IFUNC machinery results in no such function being defined, but doesn't stop include/math.h from doing the redirection resulting in the exp2m1f128 and exp10m1f128 implementations expecting to call it). Tested for x86_64 and x86, and with build-many-glibcs.py.
2024-06-17 16:31:49 +00:00
/* Return exp2(X) - 1.
Copyright (C) 2024 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
<https://www.gnu.org/licenses/>. */
#include <errno.h>
#include <fenv.h>
#include <math.h>
#include <math_private.h>
#include <math-underflow.h>
FLOAT
M_DECL_FUNC (__exp2m1) (FLOAT x)
{
if (isgreaterequal (x, M_LIT (-1.0)) && islessequal (x, M_LIT (1.0)))
{
FLOAT ret = M_SUF (__expm1) (M_MLIT (M_LN2) * x);
math_check_force_underflow (ret);
if (x != 0 && ret == 0)
__set_errno (ERANGE);
return ret;
}
else if (isgreater (x, M_MANT_DIG + M_LIT (2.0)))
{
#if defined FE_DOWNWARD || defined FE_TOWARDZERO
/* exp2m1 (MAX_EXP) should not overflow in these two rounding
modes, but exp2 does overflow. */
if (x == M_MAX_EXP)
{
int rnd_mode = fegetround ();
if (0
# ifdef FE_DOWNWARD
|| rnd_mode == FE_DOWNWARD
# endif
# ifdef FE_TOWARDZERO
|| rnd_mode == FE_TOWARDZERO
# endif
)
return M_MAX;
}
#endif
FLOAT ret = M_SUF (__ieee754_exp2) (x);
if (!isfinite (ret) && isfinite (x))
__set_errno (ERANGE);
return ret;
}
else if (isless (x, -M_MANT_DIG - M_LIT (2.0)))
return M_LIT (-1.0);
else
return M_SUF (__ieee754_exp2) (x) - M_LIT (1.0);
}
declare_mgen_alias (__exp2m1, exp2m1);