From b12a73d5bf595b7fbb73db30c6bb144078e86ef5 Mon Sep 17 00:00:00 2001
From: Szabolcs Nagy <nsz@port70.net>
Date: Tue, 11 Dec 2012 23:56:59 +0100
Subject: math: clean up inverse trigonometric functions

modifications:
* avoid unsigned->signed conversions
* removed various volatile hacks
* use FORCE_EVAL when evaluating only for side-effects
* factor out R() rational approximation instead of manual inline
* __invtrigl.h now only provides __invtrigl_R, __pio2_hi and __pio2_lo
* use 2*pio2_hi, 2*pio2_lo instead of pi_hi, pi_lo

otherwise the logic is not changed, long double versions will
need a revisit when a genaral long double cleanup happens
---
 src/math/atanl.c | 33 +++++++++++++--------------------
 1 file changed, 13 insertions(+), 20 deletions(-)

(limited to 'src/math/atanl.c')

diff --git a/src/math/atanl.c b/src/math/atanl.c
index 33ecff0d..6a480a7a 100644
--- a/src/math/atanl.c
+++ b/src/math/atanl.c
@@ -22,11 +22,6 @@ long double atanl(long double x)
 	return atan(x);
 }
 #elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
-#include "__invtrigl.h"
-
-#define ATAN_CONST      (BIAS + 65)     /* 2**65 */
-#define ATAN_LINEAR     (BIAS - 32)     /* 2**-32 */
-static const long double huge = 1.0e300;
 
 static const long double atanhi[] = {
 	 4.63647609000806116202e-01L,
@@ -81,29 +76,27 @@ long double atanl(long double x)
 	u.e = x;
 	expsign = u.xbits.expsign;
 	expt = expsign & 0x7fff;
-	if (expt >= ATAN_CONST) { /* if |x| is large, atan(x)~=pi/2 */
-		if (expt == BIAS + LDBL_MAX_EXP &&
+	if (expt >= 0x3fff + 65) { /* if |x| is large, atan(x)~=pi/2 */
+		if (expt == 0x7fff &&
 		    ((u.bits.manh&~LDBL_NBIT)|u.bits.manl)!=0)  /* NaN */
 			return x+x;
-		if (expsign > 0)
-			return  atanhi[3]+atanlo[3];
-		else
-			return -atanhi[3]-atanlo[3];
+		z = atanhi[3] + 0x1p-1000;
+		return expsign < 0 ? -z : z;
 	}
 	/* Extract the exponent and the first few bits of the mantissa. */
 	/* XXX There should be a more convenient way to do this. */
-	expman = (expt << 8) | ((u.bits.manh >> (MANH_SIZE - 9)) & 0xff);
-	if (expman < ((BIAS - 2) << 8) + 0xc0) {  /* |x| < 0.4375 */
-		if (expt < ATAN_LINEAR) {   /* if |x| is small, atanl(x)~=x */
-			/* raise inexact */
-			if (huge+x > 1.0)
-				return x;
+	expman = (expt << 8) | ((u.bits.manh >> (LDBL_MANH_SIZE - 9)) & 0xff);
+	if (expman < ((0x3fff - 2) << 8) + 0xc0) {  /* |x| < 0.4375 */
+		if (expt < 0x3fff - 32) {   /* if |x| is small, atanl(x)~=x */
+			/* raise inexact if x!=0 */
+			FORCE_EVAL(x + 0x1p1000);
+			return x;
 		}
 		id = -1;
 	} else {
 		x = fabsl(x);
-		if (expman < (BIAS << 8) + 0x30) {  /* |x| < 1.1875 */
-			if (expman < ((BIAS - 1) << 8) + 0x60) { /*  7/16 <= |x| < 11/16 */
+		if (expman < (0x3fff << 8) + 0x30) {  /* |x| < 1.1875 */
+			if (expman < ((0x3fff - 1) << 8) + 0x60) { /*  7/16 <= |x| < 11/16 */
 				id = 0;
 				x = (2.0*x-1.0)/(2.0+x);
 			} else {                                 /* 11/16 <= |x| < 19/16 */
@@ -111,7 +104,7 @@ long double atanl(long double x)
 				x = (x-1.0)/(x+1.0);
 			}
 		} else {
-			if (expman < ((BIAS + 1) << 8) + 0x38) { /* |x| < 2.4375 */
+			if (expman < ((0x3fff + 1) << 8) + 0x38) { /* |x| < 2.4375 */
 				id = 2;
 				x = (x-1.5)/(1.0+1.5*x);
 			} else {                                 /* 2.4375 <= |x| < 2^ATAN_CONST */
-- 
cgit v1.2.1