* [PATCH] Fix rounding results in lrint() & llrint() when close to 0
@ 2023-07-27 8:49 Jesse Huang
2023-07-27 9:39 ` Corinna Vinschen
0 siblings, 1 reply; 2+ messages in thread
From: Jesse Huang @ 2023-07-27 8:49 UTC (permalink / raw)
To: newlib, kito.cheng; +Cc: Jesse Huang
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
---
newlib/libm/common/s_llrint.c | 31 +++++++++++--------------------
newlib/libm/common/s_lrint.c | 31 +++++++++++--------------------
newlib/libm/common/sf_llrint.c | 6 ++----
newlib/libm/common/sf_lrint.c | 6 ++----
4 files changed, 26 insertions(+), 48 deletions(-)
diff --git a/newlib/libm/common/s_llrint.c b/newlib/libm/common/s_llrint.c
index 72452dbe9..f2c48766e 100644
--- a/newlib/libm/common/s_llrint.c
+++ b/newlib/libm/common/s_llrint.c
@@ -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)
{
diff --git a/newlib/libm/common/s_lrint.c b/newlib/libm/common/s_lrint.c
index b37f50fd4..b37b93aff 100644
--- a/newlib/libm/common/s_lrint.c
+++ b/newlib/libm/common/s_lrint.c
@@ -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)
{
diff --git a/newlib/libm/common/sf_llrint.c b/newlib/libm/common/sf_llrint.c
index 7558e89ac..905a5b21d 100644
--- a/newlib/libm/common/sf_llrint.c
+++ b/newlib/libm/common/sf_llrint.c
@@ -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
diff --git a/newlib/libm/common/sf_lrint.c b/newlib/libm/common/sf_lrint.c
index 3c58c5d10..7fe47aefb 100644
--- a/newlib/libm/common/sf_lrint.c
+++ b/newlib/libm/common/sf_lrint.c
@@ -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
--
2.31.1
^ permalink raw reply [flat|nested] 2+ messages in thread
* Re: [PATCH] Fix rounding results in lrint() & llrint() when close to 0
2023-07-27 8:49 [PATCH] Fix rounding results in lrint() & llrint() when close to 0 Jesse Huang
@ 2023-07-27 9:39 ` Corinna Vinschen
0 siblings, 0 replies; 2+ messages in thread
From: Corinna Vinschen @ 2023-07-27 9:39 UTC (permalink / raw)
To: Jesse Huang; +Cc: newlib, kito.cheng
On Jul 27 01:49, Jesse Huang via Newlib wrote:
> 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
> ---
> newlib/libm/common/s_llrint.c | 31 +++++++++++--------------------
> newlib/libm/common/s_lrint.c | 31 +++++++++++--------------------
> newlib/libm/common/sf_llrint.c | 6 ++----
> newlib/libm/common/sf_lrint.c | 6 ++----
> 4 files changed, 26 insertions(+), 48 deletions(-)
Pushed.
Thanks,
Corinna
^ permalink raw reply [flat|nested] 2+ messages in thread
end of thread, other threads:[~2023-07-27 9:39 UTC | newest]
Thread overview: 2+ messages (download: mbox.gz / follow: Atom feed)
-- links below jump to the message on this page --
2023-07-27 8:49 [PATCH] Fix rounding results in lrint() & llrint() when close to 0 Jesse Huang
2023-07-27 9:39 ` Corinna Vinschen
This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).