提交 8c44a060 编写于 作者: S Szabolcs Nagy 提交者: Rich Felker

fix scalbn when result is in the subnormal range

in nearest rounding mode scalbn could introduce double rounding error
when an intermediate value and the final result were both in the
subnormal range e.g.

  scalbn(0x1.7ffffffffffffp-1, -1073)

returned 0x1p-1073 instead of 0x1p-1074, because the intermediate
computation got rounded to 0x1.8p-1023.

with the fix an intermediate value can only be in the subnormal range
if the final result is 0 which is correct even after double rounding.
(there still can be two roundings so signals may be raised twice, but
that's only observable with trapping exceptions which is not supported.)
上级 2577b1bc
...@@ -16,11 +16,13 @@ double scalbn(double x, int n) ...@@ -16,11 +16,13 @@ double scalbn(double x, int n)
n = 1023; n = 1023;
} }
} else if (n < -1022) { } else if (n < -1022) {
y *= 0x1p-1022; /* make sure final n < -53 to avoid double
n += 1022; rounding in the subnormal range */
y *= 0x1p-1022 * 0x1p53;
n += 1022 - 53;
if (n < -1022) { if (n < -1022) {
y *= 0x1p-1022; y *= 0x1p-1022 * 0x1p53;
n += 1022; n += 1022 - 53;
if (n < -1022) if (n < -1022)
n = -1022; n = -1022;
} }
......
...@@ -16,11 +16,11 @@ float scalbnf(float x, int n) ...@@ -16,11 +16,11 @@ float scalbnf(float x, int n)
n = 127; n = 127;
} }
} else if (n < -126) { } else if (n < -126) {
y *= 0x1p-126f; y *= 0x1p-126f * 0x1p24f;
n += 126; n += 126 - 24;
if (n < -126) { if (n < -126) {
y *= 0x1p-126f; y *= 0x1p-126f * 0x1p24f;
n += 126; n += 126 - 24;
if (n < -126) if (n < -126)
n = -126; n = -126;
} }
......
...@@ -20,11 +20,11 @@ long double scalbnl(long double x, int n) ...@@ -20,11 +20,11 @@ long double scalbnl(long double x, int n)
n = 16383; n = 16383;
} }
} else if (n < -16382) { } else if (n < -16382) {
x *= 0x1p-16382L; x *= 0x1p-16382L * 0x1p113L;
n += 16382; n += 16382 - 113;
if (n < -16382) { if (n < -16382) {
x *= 0x1p-16382L; x *= 0x1p-16382L * 0x1p113L;
n += 16382; n += 16382 - 113;
if (n < -16382) if (n < -16382)
n = -16382; n = -16382;
} }
......
Markdown is supported
0% .
You are about to add 0 people to the discussion. Proceed with caution.
先完成此消息的编辑!
想要评论请 注册