This is the mail archive of the newlib@sourceware.org mailing list for the newlib project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

RE: wrong results from generic rint function (s_rint.c)


Patch attached.  Tested on i386 (i686, to be more precise).
 
The source of the problem is due to an optional step that limits how
many
bits are set in the fractional part before the primary rounding step is
done.  When the integer-part/fraction-part boundary crossed the MSB/LSB
boundary in the bit-fiddling code, the error occurred when 18 bits were
in the integral part.  The value of the fraction that could cause a
problem was basically any in which the .5 bit and any bit .25 or smaller
was also set.  (e.g. .5625, 1.75, 1.9375).
 
The apparent intent of this step of limiting the number of bits in the
fraction is to avoid double-rounding errors in implementations which
use guard bits.  (Although I can't figure out how it actually would
make any difference.  I kept it, though, in case I am missing
something.)
 
This error was simply due to a missing "else i1 = 0;".
 
However, the code also contained an inconsistency, as there are really
2 places where the number of bits in the fraction are limited.  The
first is when the int/frac boundary is near the MSB/LSB boundary, where
the error was.  It was setting 0b.001.  The second place was when it
falls only in the MSB, where it sets 0b.01.  So in addition to fixing
the error, I also changed the two fraction-limit steps to be the same,
going with the non-errored MSB-only step of 0b.01.  Not surprisingly,
sf_rint.c contained an inconsistency in this step when compared against
s_rint.c, so I have adjusted that also.  (Likely no bug in rintf(), but
it's better to have the rint() and rintf() algorithms exactly match.)
 
While fixing the problems, I also put in a few (although not complete)
comments in the code that help explain what is being done.
 
By the way, the FreeBSD code made an attempt to fix this, but I think
missed a case--with errors possibly happening with a 17-bit integer
part.
Since I'm not on the FreeBSD mailing list, if someone on this list is on
the FreeBSD list, it might be good to pass this along.
 
Craig

-----Original Message-----
From: Howland Craig D (Craig) 
Sent: Friday, March 05, 2010 10:40 AM
To: 'Ruppert'
Cc: newlib@sources.redhat.com
Subject: RE: wrong results from generic rint function (s_rint.c)

I have reproduced the problem.  I'll take a look this weekend as to
what might be wrong (assuming that nobody else jumps in and takes care
of it sooner).
Craig 

-----Original Message-----
From: newlib-owner@sourceware.org [mailto:newlib-owner@sourceware.org]
On Behalf Of Ruppert
Sent: Friday, March 05, 2010 5:33 AM
To: newlib@sources.redhat.com
Cc: dieter_ruppert@siemens.com
Subject: wrong results from generic rint function (s_rint.c)


Hi,

there seems to be a bug in the generic C implementation of rint(double) 
in s_rint.c which manifests itself in wrong rounding behaviour for
certain 
input values.

...

Regards
Dieter Ruppert

          

Attachment: rint_patch.txt
Description: rint_patch.txt

Attachment: ChangeLog.txt
Description: ChangeLog.txt


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]