Fix rounding results in lrint() & llrint() when close to 0

soft-fp should round floating pointer numbers according to the current
rounding mode. However, in the current code of lrint() and llrint(),
there are if statements before the actual rounding computation

      if(j0 < -1)
        return 0;

Where j0 is the exponent of the floating point number.

It means any number having a exponent less than -1
(i.e. interval (-0.5, 0.5)) will be rounded to 0 regardeless of the
rounding mode.

The bug already fixed in glibc in 2006 by moving the check afterwards
the rounding computation, but still persists in newlib.

This patch fixed it in a similar way to glibc
Ref Commit in glibc: 6624dbc07b5a9fb316ed188ef01f65b8eea8b47c
This commit is contained in:
Jesse Huang via Newlib 2023-07-27 01:49:41 -07:00 committed by Corinna Vinschen
parent 4fbcc8c5fe
commit 9e329b544a
4 changed files with 26 additions and 48 deletions

View File

@ -66,26 +66,17 @@ long long int
if(j0 < 20)
{
/* j0 in [-1023,19] */
if(j0 < -1)
return 0;
else
{
/* j0 in [0,19] */
/* shift amt in [0,19] */
w = TWO52[sx] + x;
t = w - TWO52[sx];
GET_HIGH_WORD(i0, t);
/* Detect the all-zeros representation of plus and
minus zero, which fails the calculation below. */
if ((i0 & ~((__int32_t)1 << 31)) == 0)
return 0;
/* After round: j0 in [0,20] */
j0 = ((i0 & 0x7ff00000) >> 20) - 1023;
i0 &= 0x000fffff;
i0 |= 0x00100000;
/* shift amt in [20,0] */
result = i0 >> (20 - j0);
}
w = TWO52[sx] + x;
t = w - TWO52[sx];
GET_HIGH_WORD(i0, t);
/* Detect the all-zeros representation of plus and
minus zero, which fails the calculation below. */
if ((i0 & ~((__int32_t)1 << 31)) == 0)
return 0;
j0 = ((i0 & 0x7ff00000) >> 20) - 1023;
i0 &= 0x000fffff;
i0 |= 0x00100000;
result = (j0 < 0 ? 0 : i0 >> (20 - j0));
}
else if (j0 < (int)(8 * sizeof (long long int)) - 1)
{

View File

@ -103,26 +103,17 @@ TWO52[2]={
if(j0 < 20)
{
/* j0 in [-1023,19] */
if(j0 < -1)
return 0;
else
{
/* j0 in [0,19] */
/* shift amt in [0,19] */
w = TWO52[sx] + x;
t = w - TWO52[sx];
GET_HIGH_WORD(i0, t);
/* Detect the all-zeros representation of plus and
minus zero, which fails the calculation below. */
if ((i0 & ~(1L << 31)) == 0)
return 0;
/* After round: j0 in [0,20] */
j0 = ((i0 & 0x7ff00000) >> 20) - 1023;
i0 &= 0x000fffff;
i0 |= 0x00100000;
/* shift amt in [20,0] */
result = i0 >> (20 - j0);
}
w = TWO52[sx] + x;
t = w - TWO52[sx];
GET_HIGH_WORD(i0, t);
/* Detect the all-zeros representation of plus and
minus zero, which fails the calculation below. */
if ((i0 & ~(1L << 31)) == 0)
return 0;
j0 = ((i0 & 0x7ff00000) >> 20) - 1023;
i0 &= 0x000fffff;
i0 |= 0x00100000;
result = (j0 < 0 ? 0 : i0 >> (20 - j0));
}
else if (j0 < (int)(8 * sizeof (long int)) - 1)
{

View File

@ -60,9 +60,7 @@ TWO23[2]={
if (j0 < (int)(sizeof (long long int) * 8) - 1)
{
if (j0 < -1)
return 0;
else if (j0 >= 23)
if (j0 >= 23)
result = (long long int) ((i0 & 0x7fffff) | 0x800000) << (j0 - 23);
else
{
@ -76,7 +74,7 @@ TWO23[2]={
j0 = ((i0 >> 23) & 0xff) - 0x7f;
i0 &= 0x7fffff;
i0 |= 0x800000;
result = i0 >> (23 - j0);
result = (j0 < 0 ? 0 : i0 >> (23 - j0));
}
}
else

View File

@ -60,9 +60,7 @@ TWO23[2]={
if (j0 < (int)(sizeof (long int) * 8) - 1)
{
if (j0 < -1)
return 0;
else if (j0 >= 23)
if (j0 >= 23)
result = (long int) ((i0 & 0x7fffff) | 0x800000) << (j0 - 23);
else
{
@ -76,7 +74,7 @@ TWO23[2]={
j0 = ((i0 >> 23) & 0xff) - 0x7f;
i0 &= 0x7fffff;
i0 |= 0x800000;
result = i0 >> (23 - j0);
result = (j0 < 0 ? 0 : i0 >> (23 - j0));
}
}
else