From: Szabolcs Nagy <nsz@port70.net>
Date: Wed, 14 Nov 2012 00:01:38 +0000 (+0100)
Subject: math: ld80 invtrig cleanups
X-Git-Tag: v0.9.8~44^2~1
X-Git-Url: https://git.librecmc.org/?a=commitdiff_plain;h=45ff9d8eb6f43b8f612d1b01ada627f5b4447c45;p=oweals%2Fmusl.git

math: ld80 invtrig cleanups
keeping only commonly used data in invtrigl
---

diff --git a/src/math/__invtrigl.c b/src/math/__invtrigl.c
index a821842c..a51330e2 100644
--- a/src/math/__invtrigl.c
+++ b/src/math/__invtrigl.c
@@ -28,9 +28,8 @@
 #include "__invtrigl.h"
 
 #if LDBL_MANT_DIG == 64 && LDBL_MAX_EXP == 16384
-/*
- * asinl() and acosl()
- */
+
+/* coefficients used in asinl() and acosl() */
 const long double
 pS0 =  1.66666666666666666631e-01L,
 pS1 = -4.16313987993683104320e-01L,
@@ -45,38 +44,9 @@ qS3 = -1.68285799854822427013e+00L,
 qS4 =  3.90699412641738801874e-01L,
 qS5 = -3.14365703596053263322e-02L;
 
-/*
- * atanl()
- */
-const long double atanhi[] = {
-	 4.63647609000806116202e-01L,
-	 7.85398163397448309628e-01L,
-	 9.82793723247329067960e-01L,
-	 1.57079632679489661926e+00L,
-};
-
-const long double atanlo[] = {
-	 1.18469937025062860669e-20L,
-	-1.25413940316708300586e-20L,
-	 2.55232234165405176172e-20L,
-	-2.50827880633416601173e-20L,
-};
-
-const long double aT[] = {
-	 3.33333333333333333017e-01L,
-	-1.99999999999999632011e-01L,
-	 1.42857142857046531280e-01L,
-	-1.11111111100562372733e-01L,
-	 9.09090902935647302252e-02L,
-	-7.69230552476207730353e-02L,
-	 6.66661718042406260546e-02L,
-	-5.88158892835030888692e-02L,
-	 5.25499891539726639379e-02L,
-	-4.70119845393155721494e-02L,
-	 4.03539201366454414072e-02L,
-	-2.91303858419364158725e-02L,
-	 1.24822046299269234080e-02L,
-};
-
+const long double pi_hi = 3.1415926535897932384626433832795L;
 const long double pi_lo = -5.01655761266833202345e-20L;
+const long double pio2_hi = 1.57079632679489661926L;
+const long double pio2_lo = -2.50827880633416601173e-20L;
+
 #endif
diff --git a/src/math/__invtrigl.h b/src/math/__invtrigl.h
index c3ad3c49..22748b68 100644
--- a/src/math/__invtrigl.h
+++ b/src/math/__invtrigl.h
@@ -32,15 +32,6 @@
 #define BIAS            (LDBL_MAX_EXP - 1)
 #define MANH_SIZE       LDBL_MANH_SIZE
 
-/* Approximation thresholds. */
-#define ASIN_LINEAR     (BIAS - 32)     /* 2**-32 */
-#define ACOS_CONST      (BIAS - 65)     /* 2**-65 */
-#define ATAN_CONST      (BIAS + 65)     /* 2**65 */
-#define ATAN_LINEAR     (BIAS - 32)     /* 2**-32 */
-
-/* 0.95 */
-#define THRESH  ((0xe666666666666666ULL>>(64-(MANH_SIZE-1)))|LDBL_NBIT)
-
 /* Constants shared by the long double inverse trig functions. */
 #define pS0     __pS0
 #define pS1     __pS1
@@ -54,56 +45,24 @@
 #define qS3     __qS3
 #define qS4     __qS4
 #define qS5     __qS5
-#define atanhi  __atanhi
-#define atanlo  __atanlo
-#define aT      __aT
+#define pi_hi   __pi_hi
 #define pi_lo   __pi_lo
+#define pio2_hi __pio2_hi
+#define pio2_lo __pio2_lo
 
-#define pio2_hi atanhi[3]
-#define pio2_lo atanlo[3]
-#define pio4_hi atanhi[1]
-
-#ifdef STRUCT_DECLS
-typedef struct longdouble {
-	uint64_t mant;
-	uint16_t expsign;
-} LONGDOUBLE;
-#else
-typedef long double LONGDOUBLE;
-#endif
-
-extern const LONGDOUBLE pS0, pS1, pS2, pS3, pS4, pS5, pS6;
-extern const LONGDOUBLE qS1, qS2, qS3, qS4, qS5;
-extern const LONGDOUBLE atanhi[], atanlo[], aT[];
-extern const LONGDOUBLE pi_lo;
-
-#ifndef STRUCT_DECLS
-static inline long double
-P(long double x)
-{
-	return (x * (pS0 + x * (pS1 + x * (pS2 + x * (pS3 + x * \
-		(pS4 + x * (pS5 + x * pS6)))))));
-}
-
-static inline long double
-Q(long double x)
-{
-	return (1.0 + x * (qS1 + x * (qS2 + x * (qS3 + x * (qS4 + x * qS5)))));
-}
+extern const long double pS0, pS1, pS2, pS3, pS4, pS5, pS6;
+extern const long double qS1, qS2, qS3, qS4, qS5;
+extern const long double pi_hi, pi_lo, pio2_hi, pio2_lo;
 
-static inline long double
-T_even(long double x)
+static long double P(long double x)
 {
-	return (aT[0] + x * (aT[2] + x * (aT[4] + x * (aT[6] + x * \
-		(aT[8] + x * (aT[10] + x * aT[12]))))));
+	return x * (pS0 + x * (pS1 + x * (pS2 + x * (pS3 +
+		x * (pS4 + x * (pS5 + x * pS6))))));
 }
 
-static inline long double
-T_odd(long double x)
+static long double Q(long double x)
 {
-	return (aT[1] + x * (aT[3] + x * (aT[5] + x * (aT[7] + x * \
-		(aT[9] + x * aT[11])))));
+	return 1.0 + x * (qS1 + x * (qS2 + x * (qS3 + x * (qS4 + x * qS5))));
 }
-#endif
 
 #endif
diff --git a/src/math/acos.c b/src/math/acos.c
index 54d266ee..0eb15bed 100644
--- a/src/math/acos.c
+++ b/src/math/acos.c
@@ -38,6 +38,7 @@
 static const double
 pi      = 3.14159265358979311600e+00, /* 0x400921FB, 0x54442D18 */
 pio2_hi = 1.57079632679489655800e+00; /* 0x3FF921FB, 0x54442D18 */
+// FIXME
 static const volatile double
 pio2_lo = 6.12323399573676603587e-17; /* 0x3C91A626, 0x33145C07 */
 static const double
diff --git a/src/math/acosl.c b/src/math/acosl.c
index cc565336..83857d49 100644
--- a/src/math/acosl.c
+++ b/src/math/acosl.c
@@ -23,9 +23,7 @@ long double acosl(long double x)
 }
 #elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
 #include "__invtrigl.h"
-
-static const long double
-pi = 3.14159265358979323846264338327950280e+00L;
+#define ACOS_CONST      (BIAS - 65)     /* 2**-65 */
 
 long double acosl(long double x)
 {
@@ -41,7 +39,8 @@ long double acosl(long double x)
 			if (expsign > 0)
 				return 0.0;  /* acos(1) = 0 */
 			else
-				return pi + 2.0 * pio2_lo;  /* acos(-1)= pi */
+				// FIXME
+				return pi_hi + 2.0 * pio2_lo;  /* acos(-1)= pi */
 		}
 		return (x - x) / (x - x);  /* acos(|x|>1) is NaN */
 	}
@@ -60,7 +59,7 @@ long double acosl(long double x)
 		s = sqrtl(z);
 		r = p / q;
 		w = r * s - pio2_lo;
-		return pi - 2.0 * (s + w);
+		return pi_hi - 2.0 * (s + w);
 	} else {                   /* x > 0.5 */
 		z = (1.0 - x) * 0.5;
 		s = sqrtl(z);
diff --git a/src/math/asinl.c b/src/math/asinl.c
index ddd807e2..7572767c 100644
--- a/src/math/asinl.c
+++ b/src/math/asinl.c
@@ -24,6 +24,10 @@ long double asinl(long double x)
 #elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
 #include "__invtrigl.h"
 static const long double huge = 1.000e+300;
+static const long double pio4_hi = 7.85398163397448309628e-01L;
+#define ASIN_LINEAR     (BIAS - 32)     /* 2**-32 */
+/* 0.95 */
+#define THRESH  ((0xe666666666666666ULL>>(64-(MANH_SIZE-1)))|LDBL_NBIT)
 
 long double asinl(long double x)
 {
diff --git a/src/math/atan2.c b/src/math/atan2.c
index 143c3834..825a38fe 100644
--- a/src/math/atan2.c
+++ b/src/math/atan2.c
@@ -39,6 +39,7 @@
 
 #include "libm.h"
 
+// FIXME
 static const volatile double
 tiny  = 1.0e-300;
 static const double
diff --git a/src/math/atan2l.c b/src/math/atan2l.c
index 0fc901c8..45cbfcc1 100644
--- a/src/math/atan2l.c
+++ b/src/math/atan2l.c
@@ -24,10 +24,8 @@ long double atan2l(long double y, long double x)
 }
 #elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
 #include "__invtrigl.h"
-static const volatile long double
-tiny = 1.0e-300;
-static const long double
-pi = 3.14159265358979323846264338327950280e+00L;
+// FIXME:
+static const volatile long double tiny = 1.0e-300;
 
 long double atan2l(long double y, long double x)
 {
@@ -55,9 +53,9 @@ long double atan2l(long double y, long double x)
 	if (expty==0 && ((uy.bits.manh&~LDBL_NBIT)|uy.bits.manl)==0) {
 		switch(m) {
 		case 0:
-		case 1: return y;        /* atan(+-0,+anything)=+-0 */
-		case 2: return  pi+tiny; /* atan(+0,-anything) = pi */
-		case 3: return -pi-tiny; /* atan(-0,-anything) =-pi */
+		case 1: return y;           /* atan(+-0,+anything)=+-0 */
+		case 2: return  pi_hi+tiny; /* atan(+0,-anything) = pi */
+		case 3: return -pi_hi-tiny; /* atan(-0,-anything) =-pi */
 		}
 	}
 	/* when x = 0 */
@@ -69,15 +67,15 @@ long double atan2l(long double y, long double x)
 			switch(m) {
 			case 0: return  pio2_hi*0.5+tiny; /* atan(+INF,+INF) */
 			case 1: return -pio2_hi*0.5-tiny; /* atan(-INF,+INF) */
-			case 2: return  1.5*pio2_hi+tiny; /*atan(+INF,-INF)*/
-			case 3: return -1.5*pio2_hi-tiny; /*atan(-INF,-INF)*/
+			case 2: return  1.5*pio2_hi+tiny; /* atan(+INF,-INF) */
+			case 3: return -1.5*pio2_hi-tiny; /* atan(-INF,-INF) */
 			}
 		} else {
 			switch(m) {
-			case 0: return  0.0;     /* atan(+...,+INF) */
-			case 1: return -0.0;     /* atan(-...,+INF) */
-			case 2: return  pi+tiny; /* atan(+...,-INF) */
-			case 3: return -pi-tiny; /* atan(-...,-INF) */
+			case 0: return  0.0;        /* atan(+...,+INF) */
+			case 1: return -0.0;        /* atan(-...,+INF) */
+			case 2: return  pi_hi+tiny; /* atan(+...,-INF) */
+			case 3: return -pi_hi-tiny; /* atan(-...,-INF) */
 			}
 		}
 	}
@@ -95,11 +93,11 @@ long double atan2l(long double y, long double x)
 	else                     /* safe to do y/x */
 		z = atanl(fabsl(y/x));
 	switch (m) {
-	case 0: return z;              /* atan(+,+) */
-	case 1: return -z;             /* atan(-,+) */
-	case 2: return pi - (z-pi_lo); /* atan(+,-) */
+	case 0: return z;               /* atan(+,+) */
+	case 1: return -z;              /* atan(-,+) */
+	case 2: return pi_hi-(z-pi_lo); /* atan(+,-) */
 	default: /* case 3 */
-		return (z-pi_lo) - pi; /* atan(-,-) */
+		return (z-pi_lo)-pi_hi; /* atan(-,-) */
 	}
 }
 #endif
diff --git a/src/math/atanl.c b/src/math/atanl.c
index 36072c17..33ecff0d 100644
--- a/src/math/atanl.c
+++ b/src/math/atanl.c
@@ -23,8 +23,53 @@ long double atanl(long double 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,
+	 7.85398163397448309628e-01L,
+	 9.82793723247329067960e-01L,
+	 1.57079632679489661926e+00L,
+};
+
+static const long double atanlo[] = {
+	 1.18469937025062860669e-20L,
+	-1.25413940316708300586e-20L,
+	 2.55232234165405176172e-20L,
+	-2.50827880633416601173e-20L,
+};
+
+static const long double aT[] = {
+	 3.33333333333333333017e-01L,
+	-1.99999999999999632011e-01L,
+	 1.42857142857046531280e-01L,
+	-1.11111111100562372733e-01L,
+	 9.09090902935647302252e-02L,
+	-7.69230552476207730353e-02L,
+	 6.66661718042406260546e-02L,
+	-5.88158892835030888692e-02L,
+	 5.25499891539726639379e-02L,
+	-4.70119845393155721494e-02L,
+	 4.03539201366454414072e-02L,
+	-2.91303858419364158725e-02L,
+	 1.24822046299269234080e-02L,
+};
+
+static long double T_even(long double x)
+{
+	return aT[0] + x * (aT[2] + x * (aT[4] + x * (aT[6] +
+		x * (aT[8] + x * (aT[10] + x * aT[12])))));
+}
+
+static long double T_odd(long double x)
+{
+	return aT[1] + x * (aT[3] + x * (aT[5] + x * (aT[7] +
+		x * (aT[9] + x * aT[11]))));
+}
+
 long double atanl(long double x)
 {
 	union IEEEl2bits u;