Attachment 229180 Details for Bug 218514 – Fix the ld128 implementations

[patch] Fix the ld128 implementations

ld128.diff (text/plain), 4.86 KB, created by Steve Kargl on 2021-10-31 22:09:41 UTC

(hide)

Description:

Filename:

MIME Type:

Creator: Steve Kargl

Created: 2021-10-31 22:09:41 UTC

Size: 4.86 KB

patch

obsolete

>diff --git a/lib/msun/ld128/s_sinpil.c b/lib/msun/ld128/s_sinpil.c
>index 39eed9b007b..cdfa2bcac3e 100644
>--- a/lib/msun/ld128/s_sinpil.c
>+++ b/lib/msun/ld128/s_sinpil.c
>@@ -1,5 +1,5 @@
> /*-
>- * Copyright (c) 2017 Steven G. Kargl
>+ * Copyright (c) 2017-2021 Steven G. Kargl
>  * All rights reserved.
>  *
>  * Redistribution and use in source and binary forms, with or without
>@@ -26,9 +26,6 @@
> 
> /*
>  * See ../src/s_sinpi.c for implementation details.
>- *
>- * FIXME:  This has not been compiled nor has it been tested for accuracy.
>- * FIXME:  This should use bit twiddling.
>  */
> 
> #include "math.h"
>@@ -68,7 +65,7 @@ sinpil(long double x)
> 			}
> 
> 			s = __kernel_sinpil(ax);
>-			return (copysignl(s, x));
>+			return (x < 0 ? -s : s);
> 		}
> 
> 		if (ax < 0.5)
>@@ -77,12 +74,18 @@ sinpil(long double x)
> 			s = __kernel_cospil(ax - 0.5);
> 		else
> 			s = __kernel_sinpil(1 - ax);
>-		return (copysignl(s, x));
>+		return (x < 0 ? -s : s);
> 	}
> 
> 	if (ax < 0x1p112) {
>-		xf = floorl(ax);
>-		ax -= xf;
>+		/* Split x = n + r with 0 <= r < 1. */
>+		xf = (ax + 0x1p112L) - 0x1p112L;        /* Integer part */
>+		ax -= xf;                               /* Remainder */
>+		if (ax < 0) {
>+			ax += 1;
>+			xf -= 1;
>+		}
>+
> 		if (ax == 0) {
> 			s = 0;
> 		} else {
>@@ -98,14 +101,14 @@ sinpil(long double x)
> 					s = __kernel_sinpil(1 - ax);
> 			}
> 
>-			if (xf > 0x1p50)
>-				xf -= 0x1p50;
>-			if (xf > 0x1p30)
>-				xf -= 0x1p30;
>+			if (xf > 0x1p64)
>+				xf -= 0x1p64;
>+			if (xf > 0x1p32)
>+				xf -= 0x1p32;
> 			ix = (uint32_t)xf;
> 			if (ix & 1) s = -s;
> 		}
>-		return (copysignl(s, x));
>+		return (x < 0 ? -s : s);
> 	}
> 
> 	if (isinf(x) || isnan(x))
>diff --git a/lib/msun/ld128/s_cospil.c b/lib/msun/ld128/s_cospil.c
>index b4bc50bb4d8..71acc4485f7 100644
>--- a/lib/msun/ld128/s_cospil.c
>+++ b/lib/msun/ld128/s_cospil.c
>@@ -1,5 +1,5 @@
> /*-
>- * Copyright (c) 2017 Steven G. Kargl
>+ * Copyright (c) 2017-2021 Steven G. Kargl
>  * All rights reserved.
>  *
>  * Redistribution and use in source and binary forms, with or without
>@@ -26,9 +26,6 @@
> 
> /*
>  * See ../src/s_cospi.c for implementation details.
>- *
>- * FIXME:  This has not been compiled nor has it been tested for accuracy.
>- * FIXME:  This should use bit twiddling.
>  */
> 
> #include "math.h"
>@@ -54,7 +51,7 @@ cospil(long double x)
> 
> 	ax = fabsl(x);
> 
>-	if (ax < 1) {
>+	if (ax <= 1) {
> 		if (ax < 0.25) {
> 			if (ax < 0x1p-60) {
> 				if ((int)x == 0)
>@@ -75,15 +72,21 @@ cospil(long double x)
> 	}
> 
> 	if (ax < 0x1p112) {
>-		xf = floorl(ax);
>-		ax -= xf;
>-		if (x < 0.5) {
>-			if (x < 0.25)
>+		/* Split x = n + r with 0 <= r < 1. */
>+		xf = (ax + 0x1p112L) - 0x1p112L;        /* Integer part */
>+		ax -= xf;                               /* Remainder */
>+		if (ax < 0) {
>+			ax += 1;
>+			xf -= 1;
>+		}
>+
>+		if (ax < 0.5) {
>+			if (ax < 0.25)
> 				c = ax == 0 ? 1 : __kernel_cospil(ax);
> 			else
> 				c = __kernel_sinpil(0.5 - ax);
> 		} else {
>-			if (x < 0.75) {
>+			if (ax < 0.75) {
> 				if (ax == 0.5)
> 					return (0);
> 				c = -__kernel_sinpil(ax - 0.5);
>@@ -91,10 +94,10 @@ cospil(long double x)
> 				c = -__kernel_cospil(1 - ax);
> 		}
> 
>-		if (xf > 0x1p50)
>-			xf -= 0x1p50;
>-		if (xf > 0x1p30)
>-			xf -= 0x1p30;
>+		if (xf > 0x1p64)
>+			xf -= 0x1p64;
>+		if (xf > 0x1p32)
>+			xf -= 0x1p32;
> 		ix = (uint32_t)xf;
> 		return (ix & 1 ? -c : c);
> 	}
>diff --git a/lib/msun/ld128/s_tanpil.c b/lib/msun/ld128/s_tanpil.c
>index 7cbbdc8a5d0..90f4aea5c62 100644
>--- a/lib/msun/ld128/s_tanpil.c
>+++ b/lib/msun/ld128/s_tanpil.c
>@@ -1,5 +1,5 @@
> /*-
>- * Copyright (c) 2017 Steven G. Kargl
>+ * Copyright (c) 2017-2021 Steven G. Kargl
>  * All rights reserved.
>  *
>  * Redistribution and use in source and binary forms, with or without
>@@ -26,9 +26,6 @@
> 
> /*
>  * See ../src/s_tanpi.c for implementation details.
>- *
>- * FIXME: This has not been compiled nor has it been tested for accuracy.
>- * FIXME: This should use bit twiddling.
>  */
> 
> #include "math.h"
>@@ -75,7 +72,7 @@ tanpil(long double x)
> 	long double ax, hi, lo, xf, t;
> 	uint32_t ix;
> 
>-	ax = fabsl(ax);
>+	ax = fabsl(x);
> 
> 	if (ax < 1) {
> 		if (ax < 0.5) {
>@@ -94,19 +91,25 @@ tanpil(long double x)
> 			return ((ax - ax) / (ax - ax));
> 		else
> 			t = -__kernel_tanpil(1 - ax);
>-		return (copysignl(t, x));
>+		return (x < 0 ? -t : t);
> 	}
> 
> 	if (ix < 0x1p112) {
>-		xf = floorl(ax);
>-		ax -= xf;
>+		/* Split x = n + r with 0 <= r < 1. */
>+		xf = (ax + 0x1p112L) - 0x1p112L;        /* Integer part */
>+		ax -= xf;                               /* Remainder */
>+		if (ax < 0) {
>+			ax += 1;
>+			xf -= 1;
>+		}
>+
> 		if (ax < 0.5)
> 			t = ax == 0 ? 0 : __kernel_tanpil(ax);
> 		else if (ax == 0.5)
> 			return ((ax - ax) / (ax - ax));
> 		else
> 			t = -__kernel_tanpil(1 - ax);
>-		return (copysignl(t, x));
>+		return (x < 0 ? -t : t);
> 	}
> 
> 	/* x = +-inf or nan. */
>@@ -114,7 +117,7 @@ tanpil(long double x)
> 		return (vzero / vzero);
> 
> 	/*
>-	 * |x| >= 0x1p53 is always an integer, so return +-0.
>+	 * |x| >= 0x1p112 is always an integer, so return +-0.
> 	 */
> 	return (copysignl(0, x));
> }

Actions: View | Diff

Attachments on bug 218514: 181622 | 182139 | 182648 | 229008 | 229010 | 229037 | 229047 | 229051 | 229108 | 229180 | 229215