Lines 63-69
Link Here
|
63 |
|
63 |
|
64 |
static __inline double pzero(double), qzero(double); |
64 |
static __inline double pzero(double), qzero(double); |
65 |
|
65 |
|
66 |
static const volatile double vone = 1, vzero = 0; |
66 |
static const volatile double vone = 1.0, vzero = 0.0; |
67 |
|
67 |
|
68 |
static const double |
68 |
static const double |
69 |
huge = 1e300, |
69 |
huge = 1e300, |
Lines 80-86
Link Here
|
80 |
S03 = 5.13546550207318111446e-07, /* 0x3EA13B54, 0xCE84D5A9 */ |
80 |
S03 = 5.13546550207318111446e-07, /* 0x3EA13B54, 0xCE84D5A9 */ |
81 |
S04 = 1.16614003333790000205e-09; /* 0x3E1408BC, 0xF4745D8F */ |
81 |
S04 = 1.16614003333790000205e-09; /* 0x3E1408BC, 0xF4745D8F */ |
82 |
|
82 |
|
83 |
static const double zero = 0.0; |
83 |
static const double zero = 0.0, qrtr = 0.25; |
84 |
|
84 |
|
85 |
double |
85 |
double |
86 |
__ieee754_j0(double x) |
86 |
__ieee754_j0(double x) |
Lines 97-103
Link Here
|
97 |
c = cos(x); |
97 |
c = cos(x); |
98 |
ss = s-c; |
98 |
ss = s-c; |
99 |
cc = s+c; |
99 |
cc = s+c; |
100 |
if(ix<0x7fe00000) { /* make sure x+x not overflow */ |
100 |
if(ix<0x7fe00000) { /* Make sure x+x does not overflow. */ |
101 |
z = -cos(x+x); |
101 |
z = -cos(x+x); |
102 |
if ((s*c)<zero) cc = z/ss; |
102 |
if ((s*c)<zero) cc = z/ss; |
103 |
else ss = z/cc; |
103 |
else ss = z/cc; |
Lines 123-131
Link Here
|
123 |
r = z*(R02+z*(R03+z*(R04+z*R05))); |
123 |
r = z*(R02+z*(R03+z*(R04+z*R05))); |
124 |
s = one+z*(S01+z*(S02+z*(S03+z*S04))); |
124 |
s = one+z*(S01+z*(S02+z*(S03+z*S04))); |
125 |
if(ix < 0x3FF00000) { /* |x| < 1.00 */ |
125 |
if(ix < 0x3FF00000) { /* |x| < 1.00 */ |
126 |
return one + z*(-0.25+(r/s)); |
126 |
return one + z*((r/s)-qrtr); |
127 |
} else { |
127 |
} else { |
128 |
u = 0.5*x; |
128 |
u = x/2; |
129 |
return((one+u)*(one-u)+z*(r/s)); |
129 |
return((one+u)*(one-u)+z*(r/s)); |
130 |
} |
130 |
} |
131 |
} |
131 |
} |
Lines 374-379
Link Here
|
374 |
static __inline double |
374 |
static __inline double |
375 |
qzero(double x) |
375 |
qzero(double x) |
376 |
{ |
376 |
{ |
|
|
377 |
static const double eighth = 0.125; |
377 |
const double *p,*q; |
378 |
const double *p,*q; |
378 |
double s,r,z; |
379 |
double s,r,z; |
379 |
int32_t ix; |
380 |
int32_t ix; |
Lines 386-390
Link Here
|
386 |
z = one/(x*x); |
387 |
z = one/(x*x); |
387 |
r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5])))); |
388 |
r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5])))); |
388 |
s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*(q[4]+z*q[5]))))); |
389 |
s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*(q[4]+z*q[5]))))); |
389 |
return (-.125 + r/s)/x; |
390 |
return (r/s-eighth)/x; |
390 |
} |
391 |
} |