提交 d8a7619e 编写于 作者: S Szabolcs Nagy

math: tgammal.c fixes

this is not a full rewrite just fixes to the special case logic:
+-0 and non-integer x<INT_MIN inputs incorrectly raised invalid
exception and for +-0 the return value was wrong

so integer test and odd/even test for negative inputs are changed
and a useless overflow test was removed
上级 e42a977f
...@@ -205,42 +205,36 @@ static long double stirf(long double x) ...@@ -205,42 +205,36 @@ static long double stirf(long double x)
long double tgammal(long double x) long double tgammal(long double x)
{ {
long double p, q, z; long double p, q, z;
int i, sign;
if (isnan(x)) if (!isfinite(x))
return NAN; return x + INFINITY;
if (x == INFINITY)
return INFINITY;
if (x == -INFINITY)
return x - x;
q = fabsl(x); q = fabsl(x);
if (q > 13.0) { if (q > 13.0) {
sign = 1;
if (q > MAXGAML)
goto goverf;
if (x < 0.0) { if (x < 0.0) {
p = floorl(q); p = floorl(q);
if (p == q)
return (x - x) / (x - x);
i = p;
if ((i & 1) == 0)
sign = -1;
z = q - p; z = q - p;
if (z > 0.5L) { if (z == 0)
p += 1.0; return 0 / z;
z = q - p; if (q > MAXGAML) {
} z = 0;
z = q * sinl(PIL * z); } else {
z = fabsl(z) * stirf(q); if (z > 0.5) {
if (z <= PIL/LDBL_MAX) { p += 1.0;
goverf: z = q - p;
return sign * INFINITY; }
z = q * sinl(PIL * z);
z = fabsl(z) * stirf(q);
z = PIL/z;
} }
z = PIL/z; if (0.5 * p == floorl(q * 0.5))
z = -z;
} else if (x > MAXGAML) {
z = x * 0x1p16383L;
} else { } else {
z = stirf(x); z = stirf(x);
} }
return sign * z; return z;
} }
z = 1.0; z = 1.0;
...@@ -268,8 +262,9 @@ goverf: ...@@ -268,8 +262,9 @@ goverf:
return z; return z;
small: small:
if (x == 0.0) /* z==1 if x was originally +-0 */
return (x - x) / (x - x); if (x == 0 && z != 1)
return x / x;
if (x < 0.0) { if (x < 0.0) {
x = -x; x = -x;
q = z / (x * __polevll(x, SN, 8)); q = z / (x * __polevll(x, SN, 8));
......
Markdown is supported
0% .
You are about to add 0 people to the discussion. Proceed with caution.
先完成此消息的编辑!
想要评论请 注册