summaryrefslogtreecommitdiff
path: root/libc/sysdeps/ieee754/flt-32/e_j1f.c
diff options
context:
space:
mode:
authorjoseph <joseph@7b3dc134-2b1b-0410-93df-9e9f96275f8d>2011-10-25 00:37:10 +0000
committerjoseph <joseph@7b3dc134-2b1b-0410-93df-9e9f96275f8d>2011-10-25 00:37:10 +0000
commit4bbe4e2185c5484328182720ff7b3bb4f9593bff (patch)
treecd67e40a74928c0f58d4f5b79d2e260e4099fee7 /libc/sysdeps/ieee754/flt-32/e_j1f.c
parent91b4be71461f78cabe1fb5f164cea71b60e9e98a (diff)
Merge changes between r15223 and r15532 from /fsf/trunk.
git-svn-id: svn://svn.eglibc.org/trunk@15545 7b3dc134-2b1b-0410-93df-9e9f96275f8d
Diffstat (limited to 'libc/sysdeps/ieee754/flt-32/e_j1f.c')
-rw-r--r--libc/sysdeps/ieee754/flt-32/e_j1f.c239
1 files changed, 65 insertions, 174 deletions
diff --git a/libc/sysdeps/ieee754/flt-32/e_j1f.c b/libc/sysdeps/ieee754/flt-32/e_j1f.c
index 71bb2515a..bb335a740 100644
--- a/libc/sysdeps/ieee754/flt-32/e_j1f.c
+++ b/libc/sysdeps/ieee754/flt-32/e_j1f.c
@@ -13,24 +13,12 @@
* ====================================================
*/
-#if defined(LIBM_SCCS) && !defined(lint)
-static char rcsid[] = "$NetBSD: e_j1f.c,v 1.4 1995/05/10 20:45:31 jtc Exp $";
-#endif
-
#include "math.h"
#include "math_private.h"
-#ifdef __STDC__
static float ponef(float), qonef(float);
-#else
-static float ponef(), qonef();
-#endif
-#ifdef __STDC__
static const float
-#else
-static float
-#endif
huge = 1e30,
one = 1.0,
invsqrtpi= 5.6418961287e-01, /* 0x3f106ebb */
@@ -46,25 +34,17 @@ s03 = 1.1771846857e-06, /* 0x359dffc2 */
s04 = 5.0463624390e-09, /* 0x31ad6446 */
s05 = 1.2354227016e-11; /* 0x2d59567e */
-#ifdef __STDC__
static const float zero = 0.0;
-#else
-static float zero = 0.0;
-#endif
-#ifdef __STDC__
- float __ieee754_j1f(float x)
-#else
- float __ieee754_j1f(x)
- float x;
-#endif
+float
+__ieee754_j1f(float x)
{
float z, s,c,ss,cc,r,u,v,y;
int32_t hx,ix;
GET_FLOAT_WORD(hx,x);
ix = hx&0x7fffffff;
- if(ix>=0x7f800000) return one/x;
+ if(__builtin_expect(ix>=0x7f800000, 0)) return one/x;
y = fabsf(x);
if(ix >= 0x40000000) { /* |x| >= 2.0 */
__sincosf (y, &s, &c);
@@ -73,7 +53,7 @@ static float zero = 0.0;
if(ix<0x7f000000) { /* make sure y+y not overflow */
z = __cosf(y+y);
if ((s*c)>zero) cc = z/ss;
- else ss = z/cc;
+ else ss = z/cc;
}
/*
* j1(x) = 1/sqrt(pi) * (P(1,x)*cc - Q(1,x)*ss) / sqrt(x)
@@ -85,9 +65,9 @@ static float zero = 0.0;
z = invsqrtpi*(u*cc-v*ss)/__ieee754_sqrtf(y);
}
if(hx<0) return -z;
- else return z;
+ else return z;
}
- if(ix<0x32000000) { /* |x|<2**-27 */
+ if(__builtin_expect(ix<0x32000000, 0)) { /* |x|<2**-27 */
if(huge+x>one) return (float)0.5*x;/* inexact if x!=0 necessary */
}
z = x*x;
@@ -96,23 +76,16 @@ static float zero = 0.0;
r *= x;
return(x*(float)0.5+r/s);
}
+strong_alias (__ieee754_j1f, __j1f_finite)
-#ifdef __STDC__
static const float U0[5] = {
-#else
-static float U0[5] = {
-#endif
-1.9605709612e-01, /* 0xbe48c331 */
5.0443872809e-02, /* 0x3d4e9e3c */
-1.9125689287e-03, /* 0xbafaaf2a */
2.3525259166e-05, /* 0x37c5581c */
-9.1909917899e-08, /* 0xb3c56003 */
};
-#ifdef __STDC__
static const float V0[5] = {
-#else
-static float V0[5] = {
-#endif
1.9916731864e-02, /* 0x3ca3286a */
2.0255257550e-04, /* 0x3954644b */
1.3560879779e-06, /* 0x35b602d4 */
@@ -120,73 +93,67 @@ static float V0[5] = {
1.6655924903e-11, /* 0x2d9281cf */
};
-#ifdef __STDC__
- float __ieee754_y1f(float x)
-#else
- float __ieee754_y1f(x)
- float x;
-#endif
+float
+__ieee754_y1f(float x)
{
float z, s,c,ss,cc,u,v;
int32_t hx,ix;
GET_FLOAT_WORD(hx,x);
- ix = 0x7fffffff&hx;
+ ix = 0x7fffffff&hx;
/* if Y1(NaN) is NaN, Y1(-inf) is NaN, Y1(inf) is 0 */
- if(ix>=0x7f800000) return one/(x+x*x);
- if(ix==0) return -HUGE_VALF+x; /* -inf and overflow exception. */
- if(hx<0) return zero/(zero*x);
- if(ix >= 0x40000000) { /* |x| >= 2.0 */
+ if(__builtin_expect(ix>=0x7f800000, 0)) return one/(x+x*x);
+ if(__builtin_expect(ix==0, 0))
+ return -HUGE_VALF+x; /* -inf and overflow exception. */
+ if(__builtin_expect(hx<0, 0)) return zero/(zero*x);
+ if(ix >= 0x40000000) { /* |x| >= 2.0 */
__sincosf (x, &s, &c);
- ss = -s-c;
- cc = s-c;
- if(ix<0x7f000000) { /* make sure x+x not overflow */
- z = __cosf(x+x);
- if ((s*c)>zero) cc = z/ss;
- else ss = z/cc;
- }
- /* y1(x) = sqrt(2/(pi*x))*(p1(x)*sin(x0)+q1(x)*cos(x0))
- * where x0 = x-3pi/4
- * Better formula:
- * cos(x0) = cos(x)cos(3pi/4)+sin(x)sin(3pi/4)
- * = 1/sqrt(2) * (sin(x) - cos(x))
- * sin(x0) = sin(x)cos(3pi/4)-cos(x)sin(3pi/4)
- * = -1/sqrt(2) * (cos(x) + sin(x))
- * To avoid cancellation, use
- * sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))
- * to compute the worse one.
- */
- if(ix>0x48000000) z = (invsqrtpi*ss)/__ieee754_sqrtf(x);
- else {
- u = ponef(x); v = qonef(x);
- z = invsqrtpi*(u*ss+v*cc)/__ieee754_sqrtf(x);
- }
- return z;
- }
- if(ix<=0x24800000) { /* x < 2**-54 */
- return(-tpi/x);
- }
- z = x*x;
- u = U0[0]+z*(U0[1]+z*(U0[2]+z*(U0[3]+z*U0[4])));
- v = one+z*(V0[0]+z*(V0[1]+z*(V0[2]+z*(V0[3]+z*V0[4]))));
- return(x*(u/v) + tpi*(__ieee754_j1f(x)*__ieee754_logf(x)-one/x));
+ ss = -s-c;
+ cc = s-c;
+ if(ix<0x7f000000) { /* make sure x+x not overflow */
+ z = __cosf(x+x);
+ if ((s*c)>zero) cc = z/ss;
+ else ss = z/cc;
+ }
+ /* y1(x) = sqrt(2/(pi*x))*(p1(x)*sin(x0)+q1(x)*cos(x0))
+ * where x0 = x-3pi/4
+ * Better formula:
+ * cos(x0) = cos(x)cos(3pi/4)+sin(x)sin(3pi/4)
+ * = 1/sqrt(2) * (sin(x) - cos(x))
+ * sin(x0) = sin(x)cos(3pi/4)-cos(x)sin(3pi/4)
+ * = -1/sqrt(2) * (cos(x) + sin(x))
+ * To avoid cancellation, use
+ * sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))
+ * to compute the worse one.
+ */
+ if(ix>0x48000000) z = (invsqrtpi*ss)/__ieee754_sqrtf(x);
+ else {
+ u = ponef(x); v = qonef(x);
+ z = invsqrtpi*(u*ss+v*cc)/__ieee754_sqrtf(x);
+ }
+ return z;
+ }
+ if(__builtin_expect(ix<=0x24800000, 0)) { /* x < 2**-54 */
+ return(-tpi/x);
+ }
+ z = x*x;
+ u = U0[0]+z*(U0[1]+z*(U0[2]+z*(U0[3]+z*U0[4])));
+ v = one+z*(V0[0]+z*(V0[1]+z*(V0[2]+z*(V0[3]+z*V0[4]))));
+ return(x*(u/v) + tpi*(__ieee754_j1f(x)*__ieee754_logf(x)-one/x));
}
+strong_alias (__ieee754_y1f, __y1f_finite)
/* For x >= 8, the asymptotic expansions of pone is
* 1 + 15/128 s^2 - 4725/2^15 s^4 - ..., where s = 1/x.
* We approximate pone by
- * pone(x) = 1 + (R/S)
+ * pone(x) = 1 + (R/S)
* where R = pr0 + pr1*s^2 + pr2*s^4 + ... + pr5*s^10
- * S = 1 + ps0*s^2 + ... + ps4*s^10
+ * S = 1 + ps0*s^2 + ... + ps4*s^10
* and
* | pone(x)-1-R/S | <= 2 ** ( -60.06)
*/
-#ifdef __STDC__
static const float pr8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
-#else
-static float pr8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
-#endif
0.0000000000e+00, /* 0x00000000 */
1.1718750000e-01, /* 0x3df00000 */
1.3239480972e+01, /* 0x4153d4ea */
@@ -194,11 +161,7 @@ static float pr8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
3.8747453613e+03, /* 0x45722bed */
7.9144794922e+03, /* 0x45f753d6 */
};
-#ifdef __STDC__
static const float ps8[5] = {
-#else
-static float ps8[5] = {
-#endif
1.1420736694e+02, /* 0x42e46a2c */
3.6509309082e+03, /* 0x45642ee5 */
3.6956207031e+04, /* 0x47105c35 */
@@ -206,11 +169,7 @@ static float ps8[5] = {
3.0804271484e+04, /* 0x46f0a88b */
};
-#ifdef __STDC__
static const float pr5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
-#else
-static float pr5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
-#endif
1.3199052094e-11, /* 0x2d68333f */
1.1718749255e-01, /* 0x3defffff */
6.8027510643e+00, /* 0x40d9b023 */
@@ -218,11 +177,7 @@ static float pr5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
5.1763616943e+02, /* 0x440168b7 */
5.2871520996e+02, /* 0x44042dc6 */
};
-#ifdef __STDC__
static const float ps5[5] = {
-#else
-static float ps5[5] = {
-#endif
5.9280597687e+01, /* 0x426d1f55 */
9.9140142822e+02, /* 0x4477d9b1 */
5.3532670898e+03, /* 0x45a74a23 */
@@ -230,11 +185,7 @@ static float ps5[5] = {
1.5040468750e+03, /* 0x44bc0180 */
};
-#ifdef __STDC__
static const float pr3[6] = {
-#else
-static float pr3[6] = {/* for x in [4.547,2.8571]=1/[0.2199,0.35001] */
-#endif
3.0250391081e-09, /* 0x314fe10d */
1.1718686670e-01, /* 0x3defffab */
3.9329774380e+00, /* 0x407bb5e7 */
@@ -242,11 +193,7 @@ static float pr3[6] = {/* for x in [4.547,2.8571]=1/[0.2199,0.35001] */
9.1055007935e+01, /* 0x42b61c2a */
4.8559066772e+01, /* 0x42423c7c */
};
-#ifdef __STDC__
static const float ps3[5] = {
-#else
-static float ps3[5] = {
-#endif
3.4791309357e+01, /* 0x420b2a4d */
3.3676245117e+02, /* 0x43a86198 */
1.0468714600e+03, /* 0x4482dbe3 */
@@ -254,11 +201,7 @@ static float ps3[5] = {
1.0378793335e+02, /* 0x42cf936c */
};
-#ifdef __STDC__
static const float pr2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
-#else
-static float pr2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
-#endif
1.0771083225e-07, /* 0x33e74ea8 */
1.1717621982e-01, /* 0x3deffa16 */
2.3685150146e+00, /* 0x401795c0 */
@@ -266,11 +209,7 @@ static float pr2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
1.7693971634e+01, /* 0x418d8d41 */
5.0735230446e+00, /* 0x40a25a4d */
};
-#ifdef __STDC__
static const float ps2[5] = {
-#else
-static float ps2[5] = {
-#endif
2.1436485291e+01, /* 0x41ab7dec */
1.2529022980e+02, /* 0x42fa9499 */
2.3227647400e+02, /* 0x436846c7 */
@@ -278,48 +217,36 @@ static float ps2[5] = {
8.3646392822e+00, /* 0x4105d590 */
};
-#ifdef __STDC__
- static float ponef(float x)
-#else
- static float ponef(x)
- float x;
-#endif
+static float
+ponef(float x)
{
-#ifdef __STDC__
const float *p,*q;
-#else
- float *p,*q;
-#endif
float z,r,s;
- int32_t ix;
+ int32_t ix;
GET_FLOAT_WORD(ix,x);
ix &= 0x7fffffff;
- if(ix>=0x41000000) {p = pr8; q= ps8;}
- else if(ix>=0x40f71c58){p = pr5; q= ps5;}
- else if(ix>=0x4036db68){p = pr3; q= ps3;}
- else if(ix>=0x40000000){p = pr2; q= ps2;}
- z = one/(x*x);
- r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5]))));
- s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*q[4]))));
- return one+ r/s;
+ if(ix>=0x41000000) {p = pr8; q= ps8;}
+ else if(ix>=0x40f71c58){p = pr5; q= ps5;}
+ else if(ix>=0x4036db68){p = pr3; q= ps3;}
+ else if(ix>=0x40000000){p = pr2; q= ps2;}
+ z = one/(x*x);
+ r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5]))));
+ s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*q[4]))));
+ return one+ r/s;
}
/* For x >= 8, the asymptotic expansions of qone is
* 3/8 s - 105/1024 s^3 - ..., where s = 1/x.
* We approximate pone by
- * qone(x) = s*(0.375 + (R/S))
+ * qone(x) = s*(0.375 + (R/S))
* where R = qr1*s^2 + qr2*s^4 + ... + qr5*s^10
- * S = 1 + qs1*s^2 + ... + qs6*s^12
+ * S = 1 + qs1*s^2 + ... + qs6*s^12
* and
* | qone(x)/s -0.375-R/S | <= 2 ** ( -61.13)
*/
-#ifdef __STDC__
static const float qr8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
-#else
-static float qr8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
-#endif
0.0000000000e+00, /* 0x00000000 */
-1.0253906250e-01, /* 0xbdd20000 */
-1.6271753311e+01, /* 0xc1822c8d */
@@ -327,11 +254,7 @@ static float qr8[6] = { /* for x in [inf, 8]=1/[0,0.125] */
-1.1849806641e+04, /* 0xc639273a */
-4.8438511719e+04, /* 0xc73d3683 */
};
-#ifdef __STDC__
static const float qs8[6] = {
-#else
-static float qs8[6] = {
-#endif
1.6139537048e+02, /* 0x43216537 */
7.8253862305e+03, /* 0x45f48b17 */
1.3387534375e+05, /* 0x4802bcd6 */
@@ -340,11 +263,7 @@ static float qs8[6] = {
-2.9449025000e+05, /* 0xc88fcb48 */
};
-#ifdef __STDC__
static const float qr5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
-#else
-static float qr5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
-#endif
-2.0897993405e-11, /* 0xadb7d219 */
-1.0253904760e-01, /* 0xbdd1fffe */
-8.0564479828e+00, /* 0xc100e736 */
@@ -352,11 +271,7 @@ static float qr5[6] = { /* for x in [8,4.5454]=1/[0.125,0.22001] */
-1.3731937256e+03, /* 0xc4aba633 */
-2.6124443359e+03, /* 0xc523471c */
};
-#ifdef __STDC__
static const float qs5[6] = {
-#else
-static float qs5[6] = {
-#endif
8.1276550293e+01, /* 0x42a28d98 */
1.9917987061e+03, /* 0x44f8f98f */
1.7468484375e+04, /* 0x468878f8 */
@@ -365,11 +280,7 @@ static float qs5[6] = {
-4.7191835938e+03, /* 0xc5937978 */
};
-#ifdef __STDC__
static const float qr3[6] = {
-#else
-static float qr3[6] = {/* for x in [4.547,2.8571]=1/[0.2199,0.35001] */
-#endif
-5.0783124372e-09, /* 0xb1ae7d4f */
-1.0253783315e-01, /* 0xbdd1ff5b */
-4.6101160049e+00, /* 0xc0938612 */
@@ -377,11 +288,7 @@ static float qr3[6] = {/* for x in [4.547,2.8571]=1/[0.2199,0.35001] */
-2.2824453735e+02, /* 0xc3643e9a */
-2.1921012878e+02, /* 0xc35b35cb */
};
-#ifdef __STDC__
static const float qs3[6] = {
-#else
-static float qs3[6] = {
-#endif
4.7665153503e+01, /* 0x423ea91e */
6.7386511230e+02, /* 0x4428775e */
3.3801528320e+03, /* 0x45534272 */
@@ -390,11 +297,7 @@ static float qs3[6] = {
-1.3520118713e+02, /* 0xc3073381 */
};
-#ifdef __STDC__
static const float qr2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
-#else
-static float qr2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
-#endif
-1.7838172539e-07, /* 0xb43f8932 */
-1.0251704603e-01, /* 0xbdd1f475 */
-2.7522056103e+00, /* 0xc0302423 */
@@ -402,11 +305,7 @@ static float qr2[6] = {/* for x in [2.8570,2]=1/[0.3499,0.5] */
-4.2325313568e+01, /* 0xc2294d1f */
-2.1371921539e+01, /* 0xc1aaf9b2 */
};
-#ifdef __STDC__
static const float qs2[6] = {
-#else
-static float qs2[6] = {
-#endif
2.9533363342e+01, /* 0x41ec4454 */
2.5298155212e+02, /* 0x437cfb47 */
7.5750280762e+02, /* 0x443d602e */
@@ -415,18 +314,10 @@ static float qs2[6] = {
-4.9594988823e+00, /* 0xc09eb437 */
};
-#ifdef __STDC__
- static float qonef(float x)
-#else
- static float qonef(x)
- float x;
-#endif
+static float
+qonef(float x)
{
-#ifdef __STDC__
const float *p,*q;
-#else
- float *p,*q;
-#endif
float s,r,z;
int32_t ix;
GET_FLOAT_WORD(ix,x);