math cleanup: use 1.0f instead of (float)1.0

This commit is contained in:
nsz
2012-03-13 20:24:23 +01:00
parent 9560b6b152
commit 8d0a6f7a1c
25 changed files with 96 additions and 96 deletions

View File

@ -24,12 +24,12 @@ static inline float __log1pf(float f)
{ {
float hfsq,s,z,R,w,t1,t2; float hfsq,s,z,R,w,t1,t2;
s = f/((float)2.0+f); s = f/(2.0f + f);
z = s*s; z = s*s;
w = z*z; w = z*z;
t1 = w*(Lg2+w*Lg4); t1 = w*(Lg2+w*Lg4);
t2 = z*(Lg1+w*Lg3); t2 = z*(Lg1+w*Lg3);
R = t2+t1; R = t2+t1;
hfsq = (float)0.5*f*f; hfsq = 0.5f * f * f;
return s*(hfsq+R); return s*(hfsq+R);
} }

View File

@ -36,8 +36,8 @@ float acosf(float x)
ix = hx & 0x7fffffff; ix = hx & 0x7fffffff;
if (ix >= 0x3f800000) { /* |x| >= 1 */ if (ix >= 0x3f800000) { /* |x| >= 1 */
if (ix == 0x3f800000) { /* |x| == 1 */ if (ix == 0x3f800000) { /* |x| == 1 */
if(hx>0) return 0.0; /* acos(1) = 0 */ if (hx > 0) return 0.0f; /* acos(1) = 0 */
return pi + (float)2.0*pio2_lo; /* acos(-1)= pi */ return pi + 2.0f*pio2_lo; /* acos(-1)= pi */
} }
return (x-x)/(x-x); /* acos(|x|>1) is NaN */ return (x-x)/(x-x); /* acos(|x|>1) is NaN */
} }
@ -50,17 +50,17 @@ float acosf(float x)
r = p/q; r = p/q;
return pio2_hi - (x - (pio2_lo-x*r)); return pio2_hi - (x - (pio2_lo-x*r));
} else if (hx < 0) { /* x < -0.5 */ } else if (hx < 0) { /* x < -0.5 */
z = (one+x)*(float)0.5; z = (one+x)*0.5f;
p = z*(pS0+z*(pS1+z*pS2)); p = z*(pS0+z*(pS1+z*pS2));
q = one+z*qS1; q = one+z*qS1;
s = sqrtf(z); s = sqrtf(z);
r = p/q; r = p/q;
w = r*s-pio2_lo; w = r*s-pio2_lo;
return pi - (float)2.0*(s+w); return pi - 2.0f*(s+w);
} else { /* x > 0.5 */ } else { /* x > 0.5 */
int32_t idf; int32_t idf;
z = (one-x)*(float)0.5; z = (one-x)*0.5f;
s = sqrtf(z); s = sqrtf(z);
df = s; df = s;
GET_FLOAT_WORD(idf,df); GET_FLOAT_WORD(idf,df);
@ -70,6 +70,6 @@ float acosf(float x)
q = one+z*qS1; q = one+z*qS1;
r = p/q; r = p/q;
w = r*s+c; w = r*s+c;
return (float)2.0*(df+w); return 2.0f*(df+w);
} }
} }

View File

@ -35,9 +35,9 @@ float acoshf(float x)
return 0.0; /* acosh(1) = 0 */ return 0.0; /* acosh(1) = 0 */
} else if (hx > 0x40000000) { /* 2**28 > x > 2 */ } else if (hx > 0x40000000) { /* 2**28 > x > 2 */
t = x*x; t = x*x;
return logf((float)2.0*x - one/(x+sqrtf(t-one))); return logf(2.0f*x - one/(x+sqrtf(t-one)));
} else { /* 1 < x < 2 */ } else { /* 1 < x < 2 */
t = x-one; t = x-one;
return log1pf(t + sqrtf((float)2.0*t+t*t)); return log1pf(t + sqrtf(2.0f*t+t*t));
} }
} }

View File

@ -52,7 +52,7 @@ float asinf(float x)
} }
/* 1 > |x| >= 0.5 */ /* 1 > |x| >= 0.5 */
w = one - fabsf(x); w = one - fabsf(x);
t = w*(float)0.5; t = w*0.5f;
p = t*(pS0+t*(pS1+t*pS2)); p = t*(pS0+t*(pS1+t*pS2));
q = one+t*qS1; q = one+t*qS1;
s = sqrt(t); s = sqrt(t);

View File

@ -38,7 +38,7 @@ float asinhf(float x)
w = logf(fabsf(x)) + ln2; w = logf(fabsf(x)) + ln2;
} else if (ix > 0x40000000) { /* 2**28 > |x| > 2.0 */ } else if (ix > 0x40000000) { /* 2**28 > |x| > 2.0 */
t = fabsf(x); t = fabsf(x);
w = logf((float)2.0*t + one/(sqrtf(x*x+one)+t)); w = logf(2.0f*t + one/(sqrtf(x*x+one)+t));
} else { /* 2.0 > |x| > 2**-28 */ } else { /* 2.0 > |x| > 2**-28 */
t = x*x; t = x*x;
w =log1pf(fabsf(x) + t/(one+sqrtf(one+t))); w =log1pf(fabsf(x) + t/(one+sqrtf(one+t)));

View File

@ -58,8 +58,8 @@ float atan2f(float y, float x)
switch (m) { switch (m) {
case 0: return pi_o_4+tiny; /* atan(+INF,+INF) */ case 0: return pi_o_4+tiny; /* atan(+INF,+INF) */
case 1: return -pi_o_4-tiny; /* atan(-INF,+INF) */ case 1: return -pi_o_4-tiny; /* atan(-INF,+INF) */
case 2: return (float)3.0*pi_o_4+tiny; /*atan(+INF,-INF)*/ case 2: return 3.0f*pi_o_4+tiny; /*atan(+INF,-INF)*/
case 3: return (float)-3.0*pi_o_4-tiny; /*atan(-INF,-INF)*/ case 3: return -3.0f*pi_o_4-tiny; /*atan(-INF,-INF)*/
} }
} else { } else {
switch (m) { switch (m) {
@ -77,7 +77,7 @@ float atan2f(float y, float x)
/* compute y/x */ /* compute y/x */
k = (iy-ix)>>23; k = (iy-ix)>>23;
if (k > 26) { /* |y/x| > 2**26 */ if (k > 26) { /* |y/x| > 2**26 */
z = pi_o_2+(float)0.5*pi_lo; z = pi_o_2 + 0.5f*pi_lo;
m &= 1; m &= 1;
} else if (k < -26 && hx < 0) /* 0 > |y|/x > -2**-26 */ } else if (k < -26 && hx < 0) /* 0 > |y|/x > -2**-26 */
z = 0.0; z = 0.0;

View File

@ -69,18 +69,18 @@ float atanf(float x)
if (ix < 0x3f980000) { /* |x| < 1.1875 */ if (ix < 0x3f980000) { /* |x| < 1.1875 */
if (ix < 0x3f300000) { /* 7/16 <= |x| < 11/16 */ if (ix < 0x3f300000) { /* 7/16 <= |x| < 11/16 */
id = 0; id = 0;
x = ((float)2.0*x-one)/((float)2.0+x); x = (2.0f*x - one)/(2.0f + x);
} else { /* 11/16 <= |x| < 19/16 */ } else { /* 11/16 <= |x| < 19/16 */
id = 1; id = 1;
x = (x-one)/(x+one); x = (x - one)/(x + one);
} }
} else { } else {
if (ix < 0x401c0000) { /* |x| < 2.4375 */ if (ix < 0x401c0000) { /* |x| < 2.4375 */
id = 2; id = 2;
x = (x-(float)1.5)/(one+(float)1.5*x); x = (x - 1.5f)/(one + 1.5f*x);
} else { /* 2.4375 <= |x| < 2**26 */ } else { /* 2.4375 <= |x| < 2**26 */
id = 3; id = 3;
x = -(float)1.0/x; x = -1.0f/x;
} }
} }
} }

View File

@ -34,9 +34,9 @@ float atanhf(float x)
SET_FLOAT_WORD(x, ix); SET_FLOAT_WORD(x, ix);
if (ix < 0x3f000000) { /* x < 0.5 */ if (ix < 0x3f000000) { /* x < 0.5 */
t = x+x; t = x+x;
t = (float)0.5*log1pf(t + t*x/(one-x)); t = 0.5f*log1pf(t + t*x/(one-x));
} else } else
t = (float)0.5*log1pf((x+x)/(one-x)); t = 0.5f*log1pf((x+x)/(one-x));
if (hx >= 0) if (hx >= 0)
return t; return t;
return -t; return -t;

View File

@ -27,7 +27,7 @@ float ceilf(float x)
if (j0 < 23) { if (j0 < 23) {
if (j0 < 0) { if (j0 < 0) {
/* raise inexact if x != 0 */ /* raise inexact if x != 0 */
if (huge+x > (float)0.0) { if (huge+x > 0.0f) {
/* return 0*sign(x) if |x|<1 */ /* return 0*sign(x) if |x|<1 */
if (i0 < 0) if (i0 < 0)
i0 = 0x80000000; i0 = 0x80000000;
@ -39,7 +39,7 @@ float ceilf(float x)
if ((i0&i) == 0) if ((i0&i) == 0)
return x; /* x is integral */ return x; /* x is integral */
/* raise inexact flag */ /* raise inexact flag */
if (huge+x > (float)0.0) { if (huge+x > 0.0f) {
if (i0 > 0) if (i0 > 0)
i0 += 0x00800000>>j0; i0 += 0x00800000>>j0;
i0 &= ~i; i0 &= ~i;

View File

@ -106,7 +106,7 @@ float erff(float x)
if (ix < 0x31800000) { /* |x| < 2**-28 */ if (ix < 0x31800000) { /* |x| < 2**-28 */
if (ix < 0x04000000) if (ix < 0x04000000)
/*avoid underflow */ /*avoid underflow */
return (float)0.125*((float)8.0*x+efx8*x); return 0.125f*(8.0f*x + efx8*x);
return x + efx*x; return x + efx*x;
} }
z = x*x; z = x*x;
@ -143,7 +143,7 @@ float erff(float x)
} }
GET_FLOAT_WORD(ix, x); GET_FLOAT_WORD(ix, x);
SET_FLOAT_WORD(z, ix&0xfffff000); SET_FLOAT_WORD(z, ix&0xfffff000);
r = expf(-z*z - (float)0.5625) * expf((z-x)*(z+x) + R/S); r = expf(-z*z - 0.5625f) * expf((z-x)*(z+x) + R/S);
if (hx >= 0) if (hx >= 0)
return one - r/x; return one - r/x;
return r/x - one; return r/x - one;
@ -206,7 +206,7 @@ float erfcf(float x)
} }
GET_FLOAT_WORD(ix, x); GET_FLOAT_WORD(ix, x);
SET_FLOAT_WORD(z, ix&0xfffff000); SET_FLOAT_WORD(z, ix&0xfffff000);
r = expf(-z*z - (float)0.5625) * expf((z-x)*(z+x) + R/S); r = expf(-z*z - 0.5625f) * expf((z-x)*(z+x) + R/S);
if (hx > 0) if (hx > 0)
return r/x; return r/x;
return two - r/x; return two - r/x;

View File

@ -85,11 +85,11 @@ float expf(float x)
SET_FLOAT_WORD(twopk, 0x3f800000+((k+100)<<23)); SET_FLOAT_WORD(twopk, 0x3f800000+((k+100)<<23));
c = x - t*(P1+t*P2); c = x - t*(P1+t*P2);
if (k == 0) if (k == 0)
return one - ((x*c)/(c-(float)2.0)-x); return one - ((x*c)/(c - 2.0f) - x);
y = one - ((lo-(x*c)/((float)2.0-c))-hi); y = one - ((lo - (x*c)/(2.0f - c)) - hi);
if (k < -125) if (k < -125)
return y*twopk*twom100; return y*twopk*twom100;
if (k == 128) if (k == 128)
return y*2.0F*0x1p127F; return y*2.0f*0x1p127f;
return y*twopk; return y*twopk;
} }

View File

@ -53,7 +53,7 @@ float expm1f(float x)
} }
if (xsb != 0) { /* x < -27*ln2 */ if (xsb != 0) { /* x < -27*ln2 */
/* raise inexact */ /* raise inexact */
if (x+tiny < (float)0.0) if (x+tiny < 0.0f)
return tiny-one; /* return -1 */ return tiny-one; /* return -1 */
} }
} }
@ -71,7 +71,7 @@ float expm1f(float x)
k = -1; k = -1;
} }
} else { } else {
k = invln2*x+((xsb==0)?(float)0.5:(float)-0.5); k = invln2*x + (xsb==0 ? 0.5f : -0.5f);
t = k; t = k;
hi = x - t*ln2_hi; /* t*ln2_hi is exact here */ hi = x - t*ln2_hi; /* t*ln2_hi is exact here */
lo = t*ln2_lo; lo = t*ln2_lo;
@ -85,27 +85,27 @@ float expm1f(float x)
k = 0; k = 0;
/* x is now in primary range */ /* x is now in primary range */
hfx = (float)0.5*x; hfx = 0.5f*x;
hxs = x*hfx; hxs = x*hfx;
r1 = one+hxs*(Q1+hxs*Q2); r1 = one+hxs*(Q1+hxs*Q2);
t = (float)3.0 - r1*hfx; t = 3.0f - r1*hfx;
e = hxs*((r1-t)/((float)6.0 - x*t)); e = hxs*((r1-t)/(6.0f - x*t));
if (k == 0) /* c is 0 */ if (k == 0) /* c is 0 */
return x - (x*e-hxs); return x - (x*e-hxs);
SET_FLOAT_WORD(twopk, 0x3f800000+(k<<23)); /* 2^k */ SET_FLOAT_WORD(twopk, 0x3f800000+(k<<23)); /* 2^k */
e = x*(e-c) - c; e = x*(e-c) - c;
e -= hxs; e -= hxs;
if (k == -1) if (k == -1)
return (float)0.5*(x-e) - (float)0.5; return 0.5f*(x-e) - 0.5f;
if (k == 1) { if (k == 1) {
if (x < (float)-0.25) if (x < -0.25f)
return -(float)2.0*(e-(x+(float)0.5)); return -2.0f*(e-(x+0.5f));
return one+(float)2.0*(x-e); return one + 2.0f*(x-e);
} }
if (k <= -2 || k > 56) { /* suffice to return exp(x)-1 */ if (k <= -2 || k > 56) { /* suffice to return exp(x)-1 */
y = one - (e - x); y = one - (e - x);
if (k == 128) if (k == 128)
y = y*2.0F*0x1p127F; y = y*2.0f*0x1p127f;
else else
y = y*twopk; y = y*twopk;
return y - one; return y - one;
@ -116,7 +116,7 @@ float expm1f(float x)
y = t - (e - x); y = t - (e - x);
y = y*twopk; y = y*twopk;
} else { } else {
SET_FLOAT_WORD(t, ((0x7f-k)<<23)); /* 2^-k */ SET_FLOAT_WORD(t, (0x7f-k)<<23); /* 2^-k */
y = x - (e + t); y = x - (e + t);
y += one; y += one;
y = y*twopk; y = y*twopk;

View File

@ -36,7 +36,7 @@ float floorf(float x)
if (j0 < 23) { if (j0 < 23) {
if (j0 < 0) { /* |x| < 1 */ if (j0 < 0) { /* |x| < 1 */
/* raise inexact if x != 0 */ /* raise inexact if x != 0 */
if (huge+x > (float)0.0) { if (huge+x > 0.0f) {
if (i0 >= 0) /* x >= 0 */ if (i0 >= 0) /* x >= 0 */
i0 = 0; i0 = 0;
else if ((i0&0x7fffffff) != 0) else if ((i0&0x7fffffff) != 0)
@ -47,7 +47,7 @@ float floorf(float x)
if ((i0&i) == 0) if ((i0&i) == 0)
return x; /* x is integral */ return x; /* x is integral */
/* raise inexact flag */ /* raise inexact flag */
if (huge+x > (float)0.0) { if (huge+x > 0.0f) {
if (i0 < 0) if (i0 < 0)
i0 += 0x00800000>>j0; i0 += 0x00800000>>j0;
i0 &= ~i; i0 &= ~i;

View File

@ -74,16 +74,16 @@ float j0f(float x)
if (huge+x > one) { if (huge+x > one) {
if (ix < 0x32000000) /* |x| < 2**-27 */ if (ix < 0x32000000) /* |x| < 2**-27 */
return one; return one;
return one - (float)0.25*x*x; return one - 0.25f*x*x;
} }
} }
z = x*x; z = x*x;
r = z*(R02+z*(R03+z*(R04+z*R05))); r = z*(R02+z*(R03+z*(R04+z*R05)));
s = one+z*(S01+z*(S02+z*(S03+z*S04))); s = one+z*(S01+z*(S02+z*(S03+z*S04)));
if(ix < 0x3F800000) { /* |x| < 1.00 */ if(ix < 0x3F800000) { /* |x| < 1.00 */
return one + z*((float)-0.25+(r/s)); return one + z*(-0.25f + (r/s));
} else { } else {
u = (float)0.5*x; u = 0.5f*x;
return (one+u)*(one-u) + z*(r/s); return (one+u)*(one-u) + z*(r/s);
} }
} }
@ -343,5 +343,5 @@ static float qzerof(float x)
z = one/(x*x); z = one/(x*x);
r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5])))); 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]+z*q[5]))))); s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*(q[4]+z*q[5])))));
return (-(float).125 + r/s)/x; return (-.125f + r/s)/x;
} }

View File

@ -75,13 +75,13 @@ float j1f(float x)
if (ix < 0x32000000) { /* |x| < 2**-27 */ if (ix < 0x32000000) { /* |x| < 2**-27 */
/* raise inexact if x!=0 */ /* raise inexact if x!=0 */
if (huge+x > one) if (huge+x > one)
return (float)0.5*x; return 0.5f*x;
} }
z = x*x; z = x*x;
r = z*(r00+z*(r01+z*(r02+z*r03))); r = z*(r00+z*(r01+z*(r02+z*r03)));
s = one+z*(s01+z*(s02+z*(s03+z*(s04+z*s05)))); s = one+z*(s01+z*(s02+z*(s03+z*(s04+z*s05))));
r *= x; r *= x;
return x*(float)0.5 + r/s; return 0.5f*x + r/s;
} }
static const float U0[5] = { static const float U0[5] = {
@ -338,5 +338,5 @@ static float qonef(float x)
z = one/(x*x); z = one/(x*x);
r = p[0]+z*(p[1]+z*(p[2]+z*(p[3]+z*(p[4]+z*p[5])))); 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]+z*q[5]))))); s = one+z*(q[0]+z*(q[1]+z*(q[2]+z*(q[3]+z*(q[4]+z*q[5])))));
return ((float).375 + r/s)/x; return (.375f + r/s)/x;
} }

View File

@ -64,7 +64,7 @@ float jnf(int n, float x)
if (n > 33) /* underflow */ if (n > 33) /* underflow */
b = zero; b = zero;
else { else {
temp = x*(float)0.5; temp = 0.5f * x;
b = temp; b = temp;
for (a=one,i=2; i<=n; i++) { for (a=one,i=2; i<=n; i++) {
a *= (float)i; /* a = n! */ a *= (float)i; /* a = n! */
@ -106,13 +106,13 @@ float jnf(int n, float x)
float q0,q1,h,tmp; float q0,q1,h,tmp;
int32_t k,m; int32_t k,m;
w = (n+n)/(float)x; w = (n+n)/x;
h = (float)2.0/(float)x; h = 2.0f/x;
z = w+h; z = w+h;
q0 = w; q0 = w;
q1 = w*z - (float)1.0; q1 = w*z - 1.0f;
k = 1; k = 1;
while (q1 < (float)1.0e9) { while (q1 < 1.0e9f) {
k += 1; k += 1;
z += h; z += h;
tmp = z*q1 - q0; tmp = z*q1 - q0;
@ -135,7 +135,7 @@ float jnf(int n, float x)
tmp = n; tmp = n;
v = two/x; v = two/x;
tmp = tmp*logf(fabsf(v*tmp)); tmp = tmp*logf(fabsf(v*tmp));
if (tmp < (float)8.8721679688e+01) { if (tmp < 88.721679688f) {
for (i=n-1,di=(float)(i+i); i>0; i--) { for (i=n-1,di=(float)(i+i); i>0; i--) {
temp = b; temp = b;
b *= di; b *= di;
@ -151,7 +151,7 @@ float jnf(int n, float x)
a = temp; a = temp;
di -= two; di -= two;
/* scale b to avoid spurious overflow */ /* scale b to avoid spurious overflow */
if (b > (float)1e10) { if (b > 1e10f) {
a /= b; a /= b;
t /= b; t /= b;
b = one; b = one;

View File

@ -104,9 +104,9 @@ static float sin_pif(float x)
*/ */
z = floorf(y); z = floorf(y);
if (z != y) { /* inexact anyway */ if (z != y) { /* inexact anyway */
y *= (float)0.5; y *= 0.5f;
y = (float)2.0*(y - floorf(y)); /* y = |x| mod 2.0 */ y = 2.0f*(y - floorf(y)); /* y = |x| mod 2.0 */
n = (int) (y*(float)4.0); n = (int)(y*4.0f);
} else { } else {
if (ix >= 0x4b800000) { if (ix >= 0x4b800000) {
y = zero; /* y must be even */ y = zero; /* y must be even */
@ -123,12 +123,12 @@ static float sin_pif(float x)
switch (n) { switch (n) {
case 0: y = __sindf(pi*y); break; case 0: y = __sindf(pi*y); break;
case 1: case 1:
case 2: y = __cosdf(pi*((float)0.5-y)); break; case 2: y = __cosdf(pi*(0.5f - y)); break;
case 3: case 3:
case 4: y = __sindf(pi*(one-y)); break; case 4: y = __sindf(pi*(one - y)); break;
case 5: case 5:
case 6: y = -__cosdf(pi*(y-(float)1.5)); break; case 6: y = -__cosdf(pi*(y - 1.5f)); break;
default: y = __sindf(pi*(y-(float)2.0)); break; default: y = __sindf(pi*(y - 2.0f)); break;
} }
return -y; return -y;
} }
@ -188,7 +188,7 @@ float lgammaf_r(float x, int *signgamp)
} else { } else {
r = zero; r = zero;
if (ix >= 0x3fdda618) { /* [1.7316,2] */ if (ix >= 0x3fdda618) { /* [1.7316,2] */
y = (float)2.0 - x; y = 2.0f - x;
i = 0; i = 0;
} else if (ix >= 0x3F9da620) { /* [1.23,1.73] */ } else if (ix >= 0x3F9da620) { /* [1.23,1.73] */
y = x - tc; y = x - tc;
@ -204,7 +204,7 @@ float lgammaf_r(float x, int *signgamp)
p1 = a0+z*(a2+z*(a4+z*(a6+z*(a8+z*a10)))); p1 = a0+z*(a2+z*(a4+z*(a6+z*(a8+z*a10))));
p2 = z*(a1+z*(a3+z*(a5+z*(a7+z*(a9+z*a11))))); p2 = z*(a1+z*(a3+z*(a5+z*(a7+z*(a9+z*a11)))));
p = y*p1+p2; p = y*p1+p2;
r += (p-(float)0.5*y); r += p - 0.5f*y;
break; break;
case 1: case 1:
z = y*y; z = y*y;
@ -218,21 +218,21 @@ float lgammaf_r(float x, int *signgamp)
case 2: case 2:
p1 = y*(u0+y*(u1+y*(u2+y*(u3+y*(u4+y*u5))))); p1 = y*(u0+y*(u1+y*(u2+y*(u3+y*(u4+y*u5)))));
p2 = one+y*(v1+y*(v2+y*(v3+y*(v4+y*v5)))); p2 = one+y*(v1+y*(v2+y*(v3+y*(v4+y*v5))));
r += (-(float)0.5*y + p1/p2); r += -0.5f*y + p1/p2;
} }
} else if (ix < 0x41000000) { /* x < 8.0 */ } else if (ix < 0x41000000) { /* x < 8.0 */
i = (int)x; i = (int)x;
y = x-(float)i; y = x - (float)i;
p = y*(s0+y*(s1+y*(s2+y*(s3+y*(s4+y*(s5+y*s6)))))); p = y*(s0+y*(s1+y*(s2+y*(s3+y*(s4+y*(s5+y*s6))))));
q = one+y*(r1+y*(r2+y*(r3+y*(r4+y*(r5+y*r6))))); q = one+y*(r1+y*(r2+y*(r3+y*(r4+y*(r5+y*r6)))));
r = half*y+p/q; r = half*y+p/q;
z = one; /* lgamma(1+s) = log(s) + lgamma(s) */ z = one; /* lgamma(1+s) = log(s) + lgamma(s) */
switch (i) { switch (i) {
case 7: z *= y + (float)6.0; /* FALLTHRU */ case 7: z *= y + 6.0f; /* FALLTHRU */
case 6: z *= y + (float)5.0; /* FALLTHRU */ case 6: z *= y + 5.0f; /* FALLTHRU */
case 5: z *= y + (float)4.0; /* FALLTHRU */ case 5: z *= y + 4.0f; /* FALLTHRU */
case 4: z *= y + (float)3.0; /* FALLTHRU */ case 4: z *= y + 3.0f; /* FALLTHRU */
case 3: z *= y + (float)2.0; /* FALLTHRU */ case 3: z *= y + 2.0f; /* FALLTHRU */
r += logf(z); r += logf(z);
break; break;
} }

View File

@ -53,8 +53,8 @@ float log10f(float x)
SET_FLOAT_WORD(x, hx|(i^0x3f800000)); /* normalize x or x/2 */ SET_FLOAT_WORD(x, hx|(i^0x3f800000)); /* normalize x or x/2 */
k += i>>23; k += i>>23;
y = (float)k; y = (float)k;
f = x - (float)1.0; f = x - 1.0f;
hfsq = (float)0.5*f*f; hfsq = 0.5f * f * f;
r = __log1pf(f); r = __log1pf(f);
// FIXME // FIXME

View File

@ -40,7 +40,7 @@ float log1pf(float x)
k = 1; k = 1;
if (hx < 0x3ed413d0) { /* 1+x < sqrt(2)+ */ if (hx < 0x3ed413d0) { /* 1+x < sqrt(2)+ */
if (ax >= 0x3f800000) { /* x <= -1.0 */ if (ax >= 0x3f800000) { /* x <= -1.0 */
if (x == (float)-1.0) if (x == -1.0f)
return -two25/zero; /* log1p(-1)=+inf */ return -two25/zero; /* log1p(-1)=+inf */
return (x-x)/(x-x); /* log1p(x<-1)=NaN */ return (x-x)/(x-x); /* log1p(x<-1)=NaN */
} }
@ -48,7 +48,7 @@ float log1pf(float x)
/* raise inexact */ /* raise inexact */
if (two25 + x > zero && ax < 0x33800000) /* |x| < 2**-24 */ if (two25 + x > zero && ax < 0x33800000) /* |x| < 2**-24 */
return x; return x;
return x - x*x*(float)0.5; return x - x*x*0.5f;
} }
if (hx > 0 || hx <= (int32_t)0xbe95f619) { /* sqrt(2)/2- <= 1+x < sqrt(2)+ */ if (hx > 0 || hx <= (int32_t)0xbe95f619) { /* sqrt(2)/2- <= 1+x < sqrt(2)+ */
k = 0; k = 0;
@ -60,11 +60,11 @@ float log1pf(float x)
return x+x; return x+x;
if (k != 0) { if (k != 0) {
if (hx < 0x5a000000) { if (hx < 0x5a000000) {
STRICT_ASSIGN(float, u, (float)1.0 + x); STRICT_ASSIGN(float, u, 1.0f + x);
GET_FLOAT_WORD(hu, u); GET_FLOAT_WORD(hu, u);
k = (hu>>23) - 127; k = (hu>>23) - 127;
/* correction term */ /* correction term */
c = k > 0 ? (float)1.0-(u-x) : x-(u-(float)1.0); c = k > 0 ? 1.0f-(u-x) : x-(u-1.0f);
c /= u; c /= u;
} else { } else {
u = x; u = x;
@ -87,9 +87,9 @@ float log1pf(float x)
SET_FLOAT_WORD(u, hu|0x3f000000); /* normalize u/2 */ SET_FLOAT_WORD(u, hu|0x3f000000); /* normalize u/2 */
hu = (0x00800000-hu)>>2; hu = (0x00800000-hu)>>2;
} }
f = u - (float)1.0; f = u - 1.0f;
} }
hfsq = (float)0.5*f*f; hfsq = 0.5f * f * f;
if (hu == 0) { /* |f| < 2**-20 */ if (hu == 0) { /* |f| < 2**-20 */
if (f == zero) { if (f == zero) {
if (k == 0) if (k == 0)
@ -97,12 +97,12 @@ float log1pf(float x)
c += k*ln2_lo; c += k*ln2_lo;
return k*ln2_hi+c; return k*ln2_hi+c;
} }
R = hfsq*((float)1.0-(float)0.66666666666666666*f); R = hfsq*(1.0f - 0.66666666666666666f * f);
if (k == 0) if (k == 0)
return f - R; return f - R;
return k*ln2_hi - ((R-(k*ln2_lo+c))-f); return k*ln2_hi - ((R-(k*ln2_lo+c))-f);
} }
s = f/((float)2.0+f); s = f/(2.0f + f);
z = s*s; z = s*s;
R = z*(Lp1+z*(Lp2+z*(Lp3+z*(Lp4+z*(Lp5+z*(Lp6+z*Lp7)))))); R = z*(Lp1+z*(Lp2+z*(Lp3+z*(Lp4+z*(Lp5+z*(Lp6+z*Lp7))))));
if (k == 0) if (k == 0)

View File

@ -51,8 +51,8 @@ float log2f(float x)
SET_FLOAT_WORD(x, hx|(i^0x3f800000)); /* normalize x or x/2 */ SET_FLOAT_WORD(x, hx|(i^0x3f800000)); /* normalize x or x/2 */
k += i>>23; k += i>>23;
y = (float)k; y = (float)k;
f = x - (float)1.0; f = x - 1.0f;
hfsq = (float)0.5*f*f; hfsq = 0.5f * f * f;
r = __log1pf(f); r = __log1pf(f);
/* /*

View File

@ -52,7 +52,7 @@ float logf(float x)
i = (ix + (0x95f64<<3)) & 0x800000; i = (ix + (0x95f64<<3)) & 0x800000;
SET_FLOAT_WORD(x, ix|(i^0x3f800000)); /* normalize x or x/2 */ SET_FLOAT_WORD(x, ix|(i^0x3f800000)); /* normalize x or x/2 */
k += i>>23; k += i>>23;
f = x - (float)1.0; f = x - 1.0f;
if ((0x007fffff & (0x8000 + ix)) < 0xc000) { /* -2**-9 <= f < 2**-9 */ if ((0x007fffff & (0x8000 + ix)) < 0xc000) { /* -2**-9 <= f < 2**-9 */
if (f == zero) { if (f == zero) {
if (k == 0) if (k == 0)
@ -60,13 +60,13 @@ float logf(float x)
dk = (float)k; dk = (float)k;
return dk*ln2_hi + dk*ln2_lo; return dk*ln2_hi + dk*ln2_lo;
} }
R = f*f*((float)0.5 - (float)0.33333333333333333*f); R = f*f*(0.5f - 0.33333333333333333f*f);
if (k == 0) if (k == 0)
return f-R; return f-R;
dk = (float)k; dk = (float)k;
return dk*ln2_hi - ((R-dk*ln2_lo)-f); return dk*ln2_hi - ((R-dk*ln2_lo)-f);
} }
s = f/((float)2.0+f); s = f/(2.0f + f);
dk = (float)k; dk = (float)k;
z = s*s; z = s*s;
i = ix-(0x6147a<<3); i = ix-(0x6147a<<3);
@ -77,7 +77,7 @@ float logf(float x)
i |= j; i |= j;
R = t2 + t1; R = t2 + t1;
if (i > 0) { if (i > 0) {
hfsq = (float)0.5*f*f; hfsq = 0.5f * f * f;
if (k == 0) if (k == 0)
return f - (hfsq-s*(hfsq+R)); return f - (hfsq-s*(hfsq+R));
return dk*ln2_hi - ((hfsq-(s*(hfsq+R)+dk*ln2_lo))-f); return dk*ln2_hi - ((hfsq-(s*(hfsq+R)+dk*ln2_lo))-f);

View File

@ -145,7 +145,7 @@ float powf(float x, float y)
/* now |1-x| is tiny <= 2**-20, suffice to compute /* now |1-x| is tiny <= 2**-20, suffice to compute
log(x) by x-x^2/2+x^3/3-x^4/4 */ log(x) by x-x^2/2+x^3/3-x^4/4 */
t = ax - 1; /* t has 20 trailing zeros */ t = ax - 1; /* t has 20 trailing zeros */
w = (t*t)*((float)0.5-t*((float)0.333333333333-t*(float)0.25)); w = (t*t)*(0.5f - t*(0.333333333333f - t*0.25f));
u = ivln2_h*t; /* ivln2_h has 16 sig. bits */ u = ivln2_h*t; /* ivln2_h has 16 sig. bits */
v = t*ivln2_l - w*ivln2; v = t*ivln2_l - w*ivln2;
t1 = u + v; t1 = u + v;
@ -193,10 +193,10 @@ float powf(float x, float y)
r = s2*s2*(L1+s2*(L2+s2*(L3+s2*(L4+s2*(L5+s2*L6))))); r = s2*s2*(L1+s2*(L2+s2*(L3+s2*(L4+s2*(L5+s2*L6)))));
r += s_l*(s_h+s); r += s_l*(s_h+s);
s2 = s_h*s_h; s2 = s_h*s_h;
t_h = (float)3.0 + s2 + r; t_h = 3.0f + s2 + r;
GET_FLOAT_WORD(is, t_h); GET_FLOAT_WORD(is, t_h);
SET_FLOAT_WORD(t_h, is & 0xfffff000); SET_FLOAT_WORD(t_h, is & 0xfffff000);
t_l = r - ((t_h - (float)3.0) - s2); t_l = r - ((t_h - 3.0f) - s2);
/* u+v = s*(1+...) */ /* u+v = s*(1+...) */
u = s_h*t_h; u = s_h*t_h;
v = s_l*t_h + t_l*s; v = s_l*t_h + t_l*s;

View File

@ -49,7 +49,7 @@ float remainderf(float x, float p)
x -= p; x -= p;
} }
} else { } else {
p_half = (float)0.5*p; p_half = 0.5f*p;
if (x > p_half) { if (x > p_half) {
x -= p; x -= p;
if (x >= p_half) if (x >= p_half)

View File

@ -19,13 +19,13 @@ float scalbf(float x, float fn)
{ {
if (isnan(x) || isnan(fn)) return x*fn; if (isnan(x) || isnan(fn)) return x*fn;
if (!isfinite(fn)) { if (!isfinite(fn)) {
if (fn > (float)0.0) if (fn > 0.0f)
return x*fn; return x*fn;
else else
return x/(-fn); return x/(-fn);
} }
if (rintf(fn) != fn) return (fn-fn)/(fn-fn); if (rintf(fn) != fn) return (fn-fn)/(fn-fn);
if ( fn > (float)65000.0) return scalbnf(x, 65000); if ( fn > 65000.0f) return scalbnf(x, 65000);
if (-fn > (float)65000.0) return scalbnf(x,-65000); if (-fn > 65000.0f) return scalbnf(x,-65000);
return scalbnf(x,(int)fn); return scalbnf(x,(int)fn);
} }

View File

@ -40,7 +40,7 @@ float sinhf(float x)
return x; return x;
t = expm1f(fabsf(x)); t = expm1f(fabsf(x));
if (ix < 0x3f800000) if (ix < 0x3f800000)
return h*((float)2.0*t - t*t/(t+one)); return h*(2.0f*t - t*t/(t+one));
return h*(t + t/(t+one)); return h*(t + t/(t+one));
} }