提交 45ff9d8e 编写于 作者: S Szabolcs Nagy

math: ld80 invtrig cleanups

keeping only commonly used data in invtrigl
上级 529a8de4
......@@ -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
......@@ -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
......@@ -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
......
......@@ -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);
......
......@@ -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)
{
......
......@@ -39,6 +39,7 @@
#include "libm.h"
// FIXME
static const volatile double
tiny = 1.0e-300;
static const double
......
......@@ -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
......@@ -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;
......
Markdown is supported
0% .
You are about to add 0 people to the discussion. Proceed with caution.
先完成此消息的编辑!
想要评论请 注册