summaryrefslogtreecommitdiff
path: root/src/math/nexttowardf.c
diff options
context:
space:
mode:
authornsz <nsz@port70.net>2012-05-06 13:08:59 +0200
committernsz <nsz@port70.net>2012-05-06 13:08:59 +0200
commit6cf865dba69bab6346dc268d9173609af36b984e (patch)
treed30e0262f1114a707bb1bc9fa03fb381e6dba324 /src/math/nexttowardf.c
parent98c9af500125df41fdb46d7e384b00982d72493a (diff)
downloadmusl-6cf865dba69bab6346dc268d9173609af36b984e.tar.gz
math: nextafter and nexttoward cleanup
make nexttoward, nexttowardf independent of long double representation. fix nextafterl: it did not raise underflow flag when the result was 0.
Diffstat (limited to 'src/math/nexttowardf.c')
-rw-r--r--src/math/nexttowardf.c81
1 files changed, 28 insertions, 53 deletions
diff --git a/src/math/nexttowardf.c b/src/math/nexttowardf.c
index c52ef3aa..821f72a5 100644
--- a/src/math/nexttowardf.c
+++ b/src/math/nexttowardf.c
@@ -1,62 +1,37 @@
-/* origin: FreeBSD /usr/src/lib/msun/src/s_nexttowardf.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.
- * ====================================================
- */
-
#include "libm.h"
-#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
-// FIXME
-#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
-#define LDBL_INFNAN_EXP (LDBL_MAX_EXP * 2 - 1)
-
float nexttowardf(float x, long double y)
{
- union IEEEl2bits uy;
- volatile float t;
- int32_t hx,ix;
-
- GET_FLOAT_WORD(hx, x);
- ix = hx & 0x7fffffff; /* |x| */
- uy.e = y;
+ union fshape ux;
+ uint32_t e;
- if (ix > 0x7f800000 ||
- (uy.bits.exp == LDBL_INFNAN_EXP &&
- ((uy.bits.manh&~LDBL_NBIT)|uy.bits.manl) != 0))
- return x + y; /* x or y is nan */
+ if (isnan(x) || isnan(y))
+ return x + y;
if (x == y)
- return (float)y; /* x=y, return y */
- if (ix == 0) { /* x == 0 */
- SET_FLOAT_WORD(x, (uy.bits.sign<<31)|1); /* return +-minsubnormal */
- /* raise underflow flag */
- t = x*x;
- if (t == x)
- return t;
- return x;
+ return y;
+ ux.value = x;
+ if (x == 0) {
+ ux.bits = 1;
+ if (signbit(y))
+ ux.bits |= 0x80000000;
+ } else if (x < y) {
+ if (signbit(x))
+ ux.bits--;
+ else
+ ux.bits++;
+ } else {
+ if (signbit(x))
+ ux.bits++;
+ else
+ ux.bits--;
}
- if (hx >= 0 ^ x < y) /* x -= ulp */
- hx--;
- else /* x += ulp */
- hx++;
- ix = hx & 0x7f800000;
- if (ix >= 0x7f800000) /* overflow */
- return x+x;
- if (ix < 0x00800000) { /* underflow */
- /* raise underflow flag */
- t = x*x;
- if (t != x) {
- SET_FLOAT_WORD(x, hx);
- return x;
- }
+ e = ux.bits & 0x7f800000;
+ /* raise overflow if ux.value is infinite and x is finite */
+ if (e == 0x7f800000)
+ return x + x;
+ /* raise underflow if ux.value is subnormal or zero */
+ if (e == 0) {
+ volatile float z = x*x + ux.value*ux.value;
}
- SET_FLOAT_WORD(x, hx);
- return x;
+ return ux.value;
}
-#endif