17 messages in org.freebsd.freebsd-standardsImplementing C99's roundf(), round(),...| From | Sent On | Attachments |
|---|---|---|
| Steve Kargl | Nov 28, 2003 4:05 pm | |
| David Schultz | Nov 29, 2003 12:13 am | |
| Steve Kargl | Nov 29, 2003 8:35 am | |
| David Schultz | Nov 30, 2003 1:42 pm | |
| Steve Kargl | Nov 30, 2003 3:53 pm | |
| Bruce Evans | Dec 1, 2003 12:05 am | |
| David Schultz | Dec 1, 2003 12:54 am | |
| Steve Kargl | Dec 1, 2003 12:35 pm | |
| Bruce Evans | Dec 1, 2003 2:57 pm | |
| Steve Kargl | Dec 1, 2003 7:39 pm | |
| Bruce Evans | Dec 1, 2003 9:57 pm | |
| David Schultz | Jan 19, 2004 1:32 pm | |
| Steve Kargl | Jan 19, 2004 2:14 pm | |
| Bruce Evans | Jan 19, 2004 8:34 pm | |
| Bruce Evans | Jan 19, 2004 10:10 pm | |
| David Schultz | Jan 21, 2004 1:18 am | |
| Bruce Evans | Jan 21, 2004 2:01 am |
| Subject: | Implementing C99's roundf(), round(), and roundl() | |
|---|---|---|
| From: | David Schultz (da...@FreeBSD.ORG) | |
| Date: | Jan 19, 2004 1:32:34 pm | |
| List: | org.freebsd.freebsd-standards | |
On Tue, Dec 02, 2003, Bruce Evans wrote:
On Mon, 1 Dec 2003, Steve Kargl wrote:
On Mon, Dec 01, 2003 at 07:05:18PM +1100, Bruce Evans wrote:
All the other corner cases need to be checked. It's possibly to check all 2^32 cases for floats (once you know the correct results).
Do you have code to do this check?
Other things to check: setting of exception flags. I'm not sure if the settings by ceil() are the right ones and the only ones.
I thought of a good way after righting the above: roundf() is supposed to be equivalent to rintf() with certain rounding, so set the rounding mode and call rintf() to determine the correct value. Unfortunately there might not be a correct rounding mode (what does round to nearest do for integers? I think it rounds to even, but roundf() requires rounding up half-way cases).
You're right on the money. The C99 standard says:
The double version of round behaves as though implemented by
#include <math.h> #include <fenv.h> #pragma STDC FENV_ACCESS ON double round(double x) { double result; fenv_t save_env; feholdexcept(&save_env); result = rint(x); if (fetestexcept(FE_INEXACT)) { fesetround(FE_TOWARDZERO); result = rint(copysign(0.5 + fabs(x), x)); } feupdateenv(&save_env); return result; }
The round functions may, but are not required to, raise the inexact exception for noninteger numeric arguments, as this implementation does.
This should be a faster implementation in that it only has one branch. The one with ceil() has a branch for the inf/NaN test, another that I added to avoid spurious inexact exceptions, and a third to test the sign of the input.
However, I think we may be able to use floor() instead of fesetround()/rint()/feupdateenv() and do something similar without changing the rounding mode:
double round(double x) { double result; fp_except_t fe;
fe = fpresetsticky(0); result = rint(x); if (fpgetsticky() & FP_X_IMP) { result = copysign(floor(fabs(x) + 0.5), x); fe |= FP_X_IMP; } fpresetsticky(fe); return (result); }
Does this seem reasonable?
By the way, the man pages refer to fpresetsticky(3), but not fpsetsticky(3). The MI ieee.h header declares only fpsetsticky(), but some architectures override these with equivalent macros for both. If you have a good idea of how this is supposed to be, it would be great if you could fix it.
Another option is to just implement the relevant parts of fenv.h. They are fairly straightforward, and there's no sense in hiding behind the brokenness of gcc's i386 backend forever. (By the way, if we did implement fenv.h, does anyone know if there are any legal issues or attribution requirements involved in using the exact round() code proposed by the C standard?)