about summary refs log tree commit diff
path: root/src/math/fma.c
diff options
context:
space:
mode:
authorSzabolcs Nagy <nsz@port70.net>2013-09-04 15:52:23 +0000
committerSzabolcs Nagy <nsz@port70.net>2013-09-05 11:30:08 +0000
commit34660d73bd0db29469d2758e1b48d2360edf3a2f (patch)
treed803e3693b46e4bfe258c4d284fc5c3f7878daae /src/math/fma.c
parent535104ab6a2d6f22098f79e7107963e3fc3448a3 (diff)
downloadmusl-34660d73bd0db29469d2758e1b48d2360edf3a2f.tar.gz
musl-34660d73bd0db29469d2758e1b48d2360edf3a2f.tar.xz
musl-34660d73bd0db29469d2758e1b48d2360edf3a2f.zip
math: fix remaining old long double code (erfl, fmal, lgammal, scalbnl)
in lgammal don't handle 1 and 2 specially, in fma use the new ldshape
union instead of ld80 one.
Diffstat (limited to 'src/math/fma.c')
-rw-r--r--src/math/fma.c38
1 files changed, 14 insertions, 24 deletions
diff --git a/src/math/fma.c b/src/math/fma.c
index 850a4a6c..1b1c1321 100644
--- a/src/math/fma.c
+++ b/src/math/fma.c
@@ -2,16 +2,6 @@
 #include "libm.h"
 
 #if LDBL_MANT_DIG==64 && LDBL_MAX_EXP==16384
-union ld80 {
-	long double x;
-	struct {
-		uint64_t m;
-		uint16_t e : 15;
-		uint16_t s : 1;
-		uint16_t pad;
-	} bits;
-};
-
 /* exact add, assumes exponent_x >= exponent_y */
 static void add(long double *hi, long double *lo, long double x, long double y)
 {
@@ -45,25 +35,25 @@ return an adjusted hi so that rounding it to double (or less) precision is corre
 */
 static long double adjust(long double hi, long double lo)
 {
-	union ld80 uhi, ulo;
+	union ldshape uhi, ulo;
 
 	if (lo == 0)
 		return hi;
-	uhi.x = hi;
-	if (uhi.bits.m & 0x3ff)
+	uhi.f = hi;
+	if (uhi.i.m & 0x3ff)
 		return hi;
-	ulo.x = lo;
-	if (uhi.bits.s == ulo.bits.s)
-		uhi.bits.m++;
+	ulo.f = lo;
+	if ((uhi.i.se & 0x8000) == (ulo.i.se & 0x8000))
+		uhi.i.m++;
 	else {
-		uhi.bits.m--;
 		/* handle underflow and take care of ld80 implicit msb */
-		if (uhi.bits.m == (uint64_t)-1/2) {
-			uhi.bits.m *= 2;
-			uhi.bits.e--;
+		if (uhi.i.m << 1 == 0) {
+			uhi.i.m = 0;
+			uhi.i.se--;
 		}
+		uhi.i.m--;
 	}
-	return uhi.x;
+	return uhi.f;
 }
 
 /* adjusted add so the result is correct when rounded to double (or less) precision */
@@ -82,9 +72,9 @@ static long double dmul(long double x, long double y)
 
 static int getexp(long double x)
 {
-	union ld80 u;
-	u.x = x;
-	return u.bits.e;
+	union ldshape u;
+	u.f = x;
+	return u.i.se & 0x7fff;
 }
 
 double fma(double x, double y, double z)